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The standard formulation of tunneling transport rests on open-boundary modeling. Most theo- 
retical studies of current-induced structural relaxations build from open-boundary calculations of 
electrostatic forces, omitting effects in the Gibbs free energy variation. A meaningful implemen- 
tation of the Born-Oppenheimer approximation (BOA) for tunneling systems requires, however, a 
consistent formulation of forces within a grand-canonical nonequilibrium (NEQ) thermodynamic 
theory, because the Friedel sum rule stipulates changes in total charge and because adiabatic trans- 
formations must be so slow that they also avoid an under-relaxed electron entropy content. In this 
paper I recast a formal NEQ thermodynamic theory for fully interacting tunneling systems, starting 
from the exact electron quantum kinetics account. I identify the operator for the Gibbs free energy 
and I define a variational NEQ thermodynamic grand potential functional which uniquely identifies 
the solution NEQ density matrix. I show that the uniqueness-of-density proof from a closely related 
Lippmann-Schwinger collision density functional theory [Phys. Rev. B 78, 165109 (2008)] makes it 
possible to express the variational NEQ thermodynamic description as a single-particle formulation 
based on universal electron-density functionals. The widely used atomistic ballistic-transport calcu- 
lations [Phys. Rev. B 52, 5335 (1995)] of the electron density are found to represent a lowest-order 
approximation to the here-presented theory of interacting tunneling. From the variational theory, I 
furthermore define thermodynamic forces which are explicitly conservative, given by a generalized 
Hellmann-Feynman theorem and different from the electrostatic forces. I demonstrate that the ther- 
modynamic forces correspond to an entropy maximization which is naturally tempered by a rigorous 
implementation of open boundary conditions. I find that these forces are, in fact, adiabatic as they 
cannot produce any thermodynamic excitation for sufficiently slow implementations of the nuclei 
dynamics. I use the variational formulation of the exact thermodynamics description to demonstrate 
that there is no inherent tendency for the electron tunneling to produce run-away nuclei-dynamics 
effects. Finally, I suggest that the here-described thermodynamic approach provides a more consis- 
tent implementation of forces, a nonequilibrium BOA for use in computational studies of interacting 
tunneling. 

PACS numbers: 73.40.Gk,71.10.-w,72.10.Bg,71.15.-m 



I. INTRODUCTION 



Understanding self-organizing functional molecular 
systems presents significant challenges for condensed- 
matter physics theory. Molecular self-assembly mech- 
anisms provide a low-cost fabrication approach* which 
mimics the molecular-recognition principle- of Nature 
and which can open technological possibilities, for ex- 
ample, atomically precise production of electronics^ and 
switchings components. The functionality is specified by 
the quantum-mechanical behavior of the electrons which 
typically respond to nonequilibrium (NEQ) conditions. 
It is natural to assume a Born-Oppenheimer approxi- 
mation (BOA) of sorts and, in practice, we solve the 
NEQ electron dynamics problem for model systems or for 
frozen nuclei coordinates within an approximation to the 
electron quantum-kinetic account (QKA)^r— Accurate 
theoretical studies have, until recently, delayed the dis- 
cussion of the NEQ forces which are exerted on the nuclei 
by transport-induced changes in the electron behavior. 
Nevertheless, it is a central result from the family of elec- 
tron density functional theory (DFT) formulations^— 
that the electron behavior is uniquely determined by 



the specific atomic configuration. Also, molecular sys- 
tems are characterized by sparseness^ i.e., have regions 
of low electron concentrations. The components are 
free to adapt their morphology to the specific environ- 
ment. To guide development of a molecular system to 
a pre-specified functionality, it is essential to develop a 
nonempiricial theory of relaxations and excitations ^2r— 
in nanoscale tunneling. 

A meaningful formulation and implementation of a 
NEQ BOA in electron QKA methods must be a corner- 
stone in the desired nonempirical theory of tunneling; the 
large ratio of nuclei to electron masses strongly suggests 
the general (but not universal) relevance of a suitable 
adiabatic behavior^ A careful analys is needed because 
the existing, equilibrium and canonical-ensemble, BOA is 
cast in terms of concepts (system staying in the evolving 
ground state (GS) and staying isolated, free of charging) 
which require a reinterpretation before use in standard 
theory formulations of NEQ tunneling. A NEQ BOA for 
a strictly electron theory of tunneling must be cast as 
a grand-canonical ensemble (GCE) thermodynamics for- 
mulation to satisfy the Gibbs free energy effects (charge 
adjustments) specified by the Friedel sum rule^ and to 
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avoid creation of thermodynamical excitations such as 
an under relaxed entropy content. In essence, the search 
for a NEQ BOA motivates a detailed discussion of the 
simplifications which must arise in an electron quantum- 
statistical mechanics (QSM) description from an assump- 
tion of adiabatic nuclei evolution. 

Recent discussions of run-away and possible water- 
wheel nuclei-dynamics effect o 48 ' — and of a speculated 
inherent nonconservative natur o 31 ' 58 ' 59 of forces accentu- 
ates the need for developing a robust framework starting 
from a consistent GCE thermodynamics theory of tun- 
neling. Sorbello has presented a gedanken-experiment ar- 
gument that forces in NEQ tunneling through a wire sys- 
tem should be nonconservativei 3 - 1 - The argument involves 
extracting an atom from the wire, moving it outside the 
wire system (taken as the Hilbert system) , and then rein- 
serting it at another point, thus seamingly avoiding the 
effect of wind forces. However, the argument is inclu- 
sive: the gedanken experiment does not consider effects 
on the dynamics produced by the extractions and inser- 
tion processes. In fact, Ref. 1571 points out that one cannot 
completely extract an atom completely from the Hilbert 
space and stay within a BOA, a complete decoupling will 
abruptly convert a resonance center into a bound state 
and therefore cause a discontinuous change in the den- 
sity of state. A complete extraction of an atom obviously 
affects the scattering properties (for example, expressed 
by phase shifts) of current-carrying electrons. Supple- 
menting and reinterpreting the discussion in Ref. 1571 . I 
note that the gedanken experiment 3 - 1 - contains an uncon- 
trolled description of the Gibbs free energy effects, for 
example, the charge adjustments that are stipulated by 
the Friedel sum rulei 5 - 6 - 

There exist clear arguments and calculations show- 
ing that modeling the dynamics of classical nuclei 
from electrostatic forces (used either as a complete 
approximatio n 3 1 > 48 ' 57 ~ ^£ or as the key starting pointS -) 
does express a nonconservative character and does lead 
not only to a possibility but also an inherent tendency for 
producing current-induced structural instabilities ! 48 ' 58 " — 
These findings discuss an effect beyond Joule heating, an 
effect which should arise also in the absence of a possi- 
ble current-induced stimulated excitation of vibrations. 
The latter effect is clearly also important, 411 for example, 
as has very recently been illustrated 5 ^ in a study of a 
current-induced limit-cycle behavior produced when fo- 
cusing the MB interaction to an electron-vibrational cou- 
pling. However, for development of an atomistic theory 
of tunneling, the most challenging findings are that the 
nonconservative electrostatic forces (defined by electron 
transport) in themselves seem capable of driving water- 
wheel 58 effects and producing a generic tendency for 
structural instabilities (even when a coupling to a phonon 
bath is also include d 60 ' 61 ). Such findings run contrary 
to the wide-spread traditional expectatio n 13 ' 28 ' 62 " — that 
NEQ tunneling can typically be analyzed as a steady- 
state transport phenomena because of the presence of an 
entropy fl ow ^MIS 



The discussions of a nonconservative nature and 
of inherent stability consequences of electrostatic 
force o 48 ' 57 " — are very important, for example, because 
such forces are now in very wide-spread usage. Use of 
electrostatic forces has become — but need not be — the 
de facto implementation choic o 40 ' 41 for nuclei relaxations 
in the widely used (yet until now, ad-hoc defined) atom- 
istic ballistic-transport codes 8 - 3 -"— for (electron) density 
calculations, here denoted BTD. However, the present- 
standard use of electrostatic forces in BTD (or in other 
infinite-system NEQ QSM) calculations is not proven 
consistent with a NEQ BOA for tunneling. 

The use of electrostatic forces in open-boundary cal- 
culations is an approximation which ignores Gibbs free 
energy interaction effects that are specified by the Friedel 
sum rule.— The use of electrostatic forces in BTD and 
in other open-boundary NEQ QSM calculations is some- 
times motivated by a referal to Ehrenfest theorem^ 7 - but 
I observe (below) that the infinite nature of an open- 
boundary description of tunneling makes an Ehrenfest- 
theorem evaluation indeterminate. Also, while there 
exists an open-boundary thermodynamical analysi o 51 ' 59 
suggesting an electrostatic nature of forces, that analysis 
contains an assumption that one can treat the number of 
left- or right-moving scattering states ('iVz,' and 'A^') as 
thermodynamical variables which are independent from 
the set nuclei positions (at fixed chemical potentials). 
This is an approximation which is not easy to reconcile 
with the Friedel sum rule 5 - 6 - that links changes in scat- 
tering potential with a net change in the integral of the 
density of state. The importance of retaining the interac- 
tion effects produced by Gibbs free energy is suggested, 
for example, by zero-temperature (equilibrium) studies of 
indirect, electron-mediated, defect interactions^"- 1 - ^ A 
change in nuclei positions, at fixed chemical potcntial(s), 
causes changes in the phase shifts, and hence concerted 
changes in the total internal energy, in the charge, and 
in the Gibbs free energy.- 1 - - 1 - 

The central role of the Gibbs free energy is well il- 
lustrated by a discussion of a thermodynamics study of 
substrate-mediated adsorbate interactions on Cu(lll). 
The Cu(lll) is a surface which supports an isotropic 
surface state and where the nonperturbative interaction 
not only can be evaluated in a GCE framewor k 101 ' 102 
but also can be compared with direct experimental 
observations j 1 - - 3 - The formal-scattering theory framework 
of Ref. llOll can easily (as also summarized in this paper) 
be used for a straightforward separate evaluation of both 
the total-internal energy and of the Gibbs free energy.^ 4 - 
However, keeping only either one of those terms leads 
to a spurious leading-order oscillatory interaction term 
(with a asymptotic decay 1/cZ,) in conflict both with the 
decay predicted in a perturbative ana lysis 1 - - - and with 
experimental observations j 1 - - 3 - Ref. llOll proceeds instead 
with a GCE evaluation (combining Mermin's equilibrium 
thermodynamics DFT— with the Harris scheme 1 - - 5 -) and 
predicts the correct I /d 2 asymptotic decay of the interac- 
tions by the subtraction of internal and Gibbs free energy 
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terms. In equilibrium this step is formally equivalent to 
pursuing a canonical-ensemble evaluation^ -100 ' 106 but 
it is not clear how one should perform a correponding 
shift to a canonical ensemble in open-boundary NEQ 
QSM or BTD studies. Instead the experience from 
the equilibrium study of surface-state mediated defect 
interactions 103 suggests an importance of retaining the 
Gibbs free energy term in an open-boundary thermody- 
namics determination of forces. 

I stress that many-body (MB) interactions in NEQ 
tunneling are known to cause changes in the statisti- 
cal distribution among MB eigenstates^&2£u & 104 ' 107 
in tunneling. The magnitude and complexity of such 
changes significantly exceed those arising in equilibrium 
cases because the Fermi-level concept absolves and gen- 
uine inelastic processes can proceed] 71 ' 77 ' 82 ' 87 The Friedel 
sum rule is, of course, also of central importance in 
correct calculations of such MB interactions In ad- 
dition, it is in NEQ tunneling a much more complex 
task to determine what possible MB states represent the 
lower- or higher-lying MB configurations of the emerging 
distributioni^ A full, exact GCE thermodynamical ac- 
count is essential to decide which such MB configurations 
must be excluded in an implementation of an adiabatic 
transformation of nuclei positions. While use of electro- 
static forces does ensure a disentanglement of nuclei and 
electron dynamics, there is no argument that it will also 
avoid creating such MB thermodynamical excitations of 
the ensemble. It is motivated to deepen the discussion of 
an adiabatic nature in NEQ tunneling. 

This paper recasts a rigorous NEQ QSM formulation of 
GCE thermodynamics of general NEQ tunneling, thereby 
permitting a precise discussion of (1) the role of the Gibbs 
free energy, (2) of a formal single-particle (SP) formu- 
lation (to which the BTD method for electron-density 
calculations is found to represent a lowest-order yet con- 
sistent approximation) and (3) of the adiabatic nature of 
relaxations in NEQ tunneling. The paper delivers a suit- 
able definition of thermodynamic forces and allows an 
implementation of a NEQ BOA for the electron QKA. 
To accomplish that, the paper expresses a formally exact 
NEQ thermodynamic theory within the electron QKA 
which, in turn, is expected to contain all boundary, en- 
tropy, and entropy-flow effects in the electron system. 
The paper proves a generalization of the variational prin- 
ciple for the NEQ thermodynamic GCE; the proof follows 
by a straightforward extension of Mermin's analysis in his 
development of the equilibrium-thermodynamic DFTJ^ 
The variational property reflects a maximization of en- 
tropy subject to the rigorous boundary conditions set by 
exact QKA determinations of the internal energy and 
Gibbs-free energy contents. The here-derived thermo- 
dynamic forces are conservative and cause the electron- 
nuclei system to reach a further entropy maximization 
(again tempered by evaluations of relevant thermody- 
namic properties). 

I note that the exact NEQ GCE thermodynamic theory 
does not provide (nor need) a uniqueness-of-density argu- 



ment (as is possible in the equilibrium CclSC ). However, 
I also note that the present NEQ thermodynamic account 
is developed in the same continuum Caroli partition 
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as was previously used in the formulation 2 * of 



a MB Lippmann-Schwinge r 109 ' 110 (LS) collision density 
functional theory (LSC-DFT); the appendix of Ref. [28| 
provides an argument that one can treat both the NEQ 
thermodynamic quantities and associated NEQ density 
matrix solutions as electron-density functionals. The SP 
formulation (and the expression of BTD as an approxi- 
mation to the formally exact GCE thermodynamics the- 
ory) builds on the uniqueness of density theorem pre- 
sented in Ref. |H 

The paper is organized as follows. The next section 
presents the background, summarizing the recent discus- 
sions of a possible nonconservative nature of forces to 
which I am here adding a perspective on the role of the 
Gibbs free energy and on adiabacity. Section III con- 
tains an overview of the formal results for thermody- 
namic quantities and a description of the nature of the 
here-derived thermodynamic forces. Section IV presents 
the partition scheme, including the Hamiltonian and the 
thermodynamic boundary conditions which are used in 
the formal collision description of NEQ interaction tun- 
neling and upon which the thermodynamic description 
rests. Section V provides a presentation of the QKA and 
formal NEQ QSM solutions, while Sec. VI presents the 
variational nature of the NEQ thermodynamic theory. 
The linking with the LSC-DFT 2 ^ and the SP formula- 
tion of the interacting thermodynamical theory is pre- 
sented in Sec. VII. Section VIII investigates the nature 
of steady-state forces and Sec. IX contains a summary 
and acknowledgment. Finally, the paper contains two 
appendices. Appendix A provides details of the proof 
for a Mermin-type variational principle, i.e., a general- 
ized Hellmann-Feynman (GHF) theorem. Appendix B 
contains a mathematical argument for formal compati- 
bility of LSC-DFT and the here-presented NEQ thermo- 
dynamic theory. 



II. BACKGROUND: FORCES, RELAXATIONS, 
AND STABILITY IN TUNNELING 



The highly successful standard DFT^ii^ i.e., GS- 
DFT, illustrates the level of detail and predictive power 
that we desire from a new nonempirical NEQ com- 
putational theory. GS-DFT rests on an equilibrium 
total-energy variational principle and regularly deliv- 
ers accurate characterizations of structure for dense 
(hard) matter problems using the generalized gradient 
approximation*^ In GS-DFT, the forces F^' GS on the 
nuclei are entirely electrostatic (hence identified by a su- 
perscript 'es') in nature and they are given by the tradi- 
tional, GS Hellmann-Feynman (HF) theorem, 
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connecting the total- internal energy C/qs = (H)gs di- 
rectly with the equilibrium forces. 

With accurate and transferable functional approxima- 
tions, the GS-DFT is ideally suited to provide a nonem- 
pirical description of structural relaxations through an 
equilibrium BOA; calculations of the density provide a di- 
rect determination of the Coulombic force from Eq. (JlJ . 
The flexibility for relaxations and the presence of im- 
portant low density regions does imply that even GS- 
DFT faces a significant challenge when trying to describe 
molecular tunneling system already in equilibrium. Nev- 
ertheless, GS-DFT has recently begun to deliver results 
for the even broader class of sparse-matter systems i 29 ' 112 
The presence of low electron density regions in such sys- 
tems implies that the dispersive and van der Waals forces 
play an important role in determining the resulting struc- 
ture, cohesion, and functionality. There are develop- 
ments both in semi-empirical^ 1 ^ and new nonlocal den- 
sity functionat^~-ii£ descriptions with accounts of dis- 
persive or van der Waals interactions. It has become pos- 
sible to obtain, for example, nonempirical descriptions of 
interactions and resulting morphology of self-organizing 
organic-molecular overlayer systems 1 -^-^ in equilibrium 
(and as described in a canonical-ensemble theory). 

However, in spite of the GS-DFT versatility, there is 
no meaningful application to working tunneling systems 
since the NEQ nature invalidates the foundation for the 
total-energy variational principle. A breakdown of the 
GS HF theorem is broadly expected. The source and 
extent of the deviation (of the forces from the derivative 
of the total internal energy) are explicitly expressed in 
this paper. 

I note that a generalized HF-type variational prop- 
erty, i.e., a GHF theorem defined by the thermody- 
namic grand potential^ permits a corresponding sim- 
plification of forces in systems described as an equilib- 
rium GCE theory. However, the GCE thermodynamics 
forces, F^ CE ^ F^' E , must also reflect a variation 
of the Gibbs free energ y 101 ' 104 The GCE thermodynam- 
ics forces can not seen being strictly electrostatic in na- 
ture when we assume fixed chemical potential(s). This 
is true even though the values (i.e., the force vectors) 
must, of course, coincide with those obtained i n a corre- 
sponding canonical-ensemble evaluation, Refs. ll0lUl06l 
A documented limitations 1 ^ of electrostatic forces for a 
GCE equilibrium (zero-temperature) case is briefly sum- 
marized in this paper. 

The BTD method 8 ^— is a widely used approach for 
atomistic calculations of the electron density (and trans- 
port) in NEQ tunneling systems described in the stan- 
dard open-boundary framework. The BTD method rests 
on a self-consistency loop in which a Landauer-type 6 ^ 
ballistic transport picture is used to calculate the electron 
density from SP scattering states and in which that den- 
sity is used to define the effective SP scattering potential 
(using the GS-DFT density functionals to describe some 
MB interaction effects). The BTD method is robust, as 
it codes SP LS solutions of the SP GCE thermodynamics 



problem and thus (automatically) satisfies conservation 
laws relevant for the effective SP potential. In particular, 
the BTD calculations of the electron density do reflect 
charge adjustments (Gibbs-free energy effects) mandated 
by the Friedel sum rule for the effective potential. For 
investigations of current-induced relaxations, the BTD 
method is typically coupled with an evaluation and im- 
plementations of electrostatic forces ; 40 ' 41 but as discussed 
above (and below), there exists no conclusive argument 
for that choice. The BTD method for electron-density 
calculations is generally seen as extremely succesful in 
delivering efficient, parameter- free calculations of tunnel- 
ing at an ab initio level of detail. Existing now in very 
efficient codings and allowing even for inclusion of disper- 
sive of vdW interactions, it is emerging as a powerful tool 
for detailed atomistic-theory studies of also molecular- 
electronics systems. I note that it is the emphasis on 
coding conservation laws, i.e., physics principles, which 
ensures the tremendious success and progress of GS-DFT 
for equilibrium systems. Since the BTD method delivers 
an impressive level of details in predictions and charac- 
terizations, there is an expectation that the BTD method 
derives from a robust, albeit (until now) not identified, 
underlying physics principle. This paper finds that the 
BTD method for electron density calculations (but not 
a coupling of BTD with electrostatic forces) can be seen 
as a lowest-order approximation for a formally exact SP 
calculation scheme for a full determination of the CGE 
thermodynamics of tunneling. 

The QKASr— provides, in principle, a full account 
of NEQ QSM but actual implementations focus on the 
electron dynamics and has (until now) in practice been 
forced to rely on empirical input for descriptions of ac- 
tual systems. The formal electron QKA gives an ex- 
act treatment of NEQ tunneling in a general external 
(configuration-specific) potential for the standard open- 
boundary modeling framework. Actual computational 
implementations of QKA methods are (of course) ap- 
proximative but, worse, they are also empirical in the 
sense that they rely on assumptions of a frozen external 
potential (i.e., an assumption of frozen nuclei positions) 
when used to describe actual systems. It is because of 
this frozen-nuclei computation approach, that we must 
now (separately) develop a suitable implementation of 
a NEQ BOA to correctly reflect the expected adiabatic 
behavior. 

The electron QKA is very powerful and useful. Of cen- 
tral importance for this paper, I note that the set of QKA 
methods includes time-evolution formulation s] 120 ' 121 and 
reformulations (here termed ReNQSM) of the NEQ QSM 
for steady-state tunnelin g 13 ' 67 ' 68 ' 92 ' 107 and for tunnel- 
ing systems with a period or pulsed time dependence 
of the external field . 122 ' 123 I note that the ReNQSM 
framework can be adapted for both time-averaged^ 3 - and 
variational 9 — determinations of interacting transport in 
nanoscale systems. The electron QKA also includes 
Wigner-distribution calculations, for example, as de- 
scribed in an early and clear discussion^ 6 - of the role of 
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open boundaries and of a need for a thermodynamical 
treatment of NEQ tunneling. 

The electron QKA has found its most common expres- 
sion in the NEQ Green function (NEGF) techniques- 
based on a diagrammatic approach and with an 
important emphasis on conserving approximations. 7 
The NEGF approach has brought both analysis and 
methods^^i^^^^^^il-^^^^^ Part of this 
formal development involves a coupling of an elec- 
tron QKA with both traditional (finite-system) TD- 
DFT^ 7 . and adaptations of TD-DFT results to in- 
finite tunneling systemsi 25 i 28 ' 126 In the context of 
nanoscale tunneling transport, there exist, for ex- 
ample, formal generalization o 65 ' 71 ^ 72 of the Landauer- 
Biittikcr formula 64 to MB interaction cases but also 
warning s 67 i 71 ' 77 i 82 ' 87 i 126 i 127 of a breakdown for an over- 
simplified SP interpretation of this formal result. The 
NEGF methods have in turn given us a detailed charac- 
terization of NEQ interaction effects in nanoscale tunnel- 
ing trancport . 32 - 39,45,47 - 49 - 50 - 66 - 70 - 71 - 75 - 79 - 81 - 82 - 88-90 - 128-133 

At the same time, I note that there is a lack of a nonem- 
pirical quality in present-day electron-QKA studies which 
do not implement effects of NEQ forces and relaxations 
(except in model studies). It is exciting that recent 
NEQ QSM calculations, using the self-consistent GW 



approximation, 



130,132,133 



have begun to free the NEGF 



method of rigid, frozen parameterizations of the inter- 
atomic coupling strengths (for actual systems). Still, 
such (electron) calculations are not yet employed to also 
predict structural relaxations which results, for example, 
when we seek to characterize the multiple (meta-)stable 
configurations in current-induced molecular switching. 4 
The absence of a nonempirical character arises not only 
from computational limitations but also from a funda- 
mental problem: there is no universally accepted ap- 
proach to implement the BOA for systems that are out 
of equilibrium. A complete QKA, treating both elec- 
trons and nuclei as quantum-mechanical objects, will ob- 
viously provide an exact determination of wind forces and 
electromigratio n 30 ! 31 ' 52 
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However, the key issue is to ensure a suitable adiabatic 
behavior^ i.e. a meaningfull NEQ BOA, when using 
electron-QKA calculations to approximately predict the 
actual relaxation. The challenge of trying to thus ad- 
dress, for example, structural relaxations and the elec- 
tromigration in an electron GCE thermodynamics deter- 
mination of adiabatic forces 
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has stimulated a sequence of method-development, anal- 
ysis, and model studies 48,52-54,57-59 

Almost all existing descriptions of tunneling forces 
starts with the definition 
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(where the subscript on the right-hand side emphasize 
that the electrostatic force is evaluated in a NEQ GCE 
description as mandated by the standard open-boundary 
formulation) and the assumption 
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I note that while the assumption of electrostatic forces 
has recently been supplemented by terms describing the 
coupling to a phonon bath, it is the electron behavior 
expressed by the Coulombic forces (j4|) which are typically 
seen (in open-boundary investigations) as the source for 
a possible structural instability in tunneling systems i 60 ! 61 

The assumption of electrostatic forces (J5J) enters, for 
example, as expressed by a choice of Lagrangian^ or 
when motivation is sought from Ehrenfest theorem^ 7 , fol- 
lowed by an Ehrenfest approximation (i.e., a separa- 
tion of slow nuclei and fast electron dynamics 5 - 2 - - — ) al- 
though the standard open-boundary (and hence infinite- 
system) modeling renders an Ehrenfest evaluation inde- 
terminate (below). Ehrenfest theorem is, of course, appli- 
cable for calculations in traditional, finite-extention time- 
dependent DFT (TD-DFT) 134 ' 135 and when the NEQ 
transport problem is recast as a finite canonical-ensemble 
problem where TD-DFT describe a driving force bal- 
anced by a coupling to a phonon heat bathi 26 i 27 How- 
ever, this paper focus on the standard open-boundary 
modeling of tunneling and will make no further discus- 
sion of the canonical-ensemble approach. All discussion 
and derivations below (sometimes implicitly) assumes the 
open-boundary modeling framework where it is impera- 
tive to use a GCE thermodynamics description. 

Electrostatic forces are a difficult assumption in a GCE 
description. A question mark on the validity of the rela- 
tion 
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must be clearly expressed in any (NEQ or equilibrium) 
GCE evaluations of adiabatic forces. The question marks 
in ([3]) and ([5]) applies in particular to open-boundary 
model framework for tunneling for a number of rea- 
sons involving mathematical rigor, the Friedel sum rule, 
and the thermodynamics which is implied by the exact 
steady-state ReNQSM^ 7 . 

First, the assumption ([5]) lacks a formal derivation for 
open, infinite systems. Like Refs. l53lJ54j|57l . I observe that 
the force description ([5]) can only be derived subject to 
an additional assumption of an absence of a velocity- 
dependent force term. This is a relevant reservation 
for one must expect that friction effects can arise not 
only through coupling to vibrations but also within an 
electron QKA, as a result of a combination of electron- 
electron scatte ring an d of the finite particle and entropy 
flow. Refs. I52il57ll59l correctly point out that there is a 
need for corrections to the use of Ehrenfest theorem^ 7 - 
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Here the subscripts on both sides emphasize that the 
proof is only available for for finite systems^ or in 
a standard GS-DFT formulation; 1 ^ 6 . Reservations are 
formulate d 52 ' 57 ' 59 when we have only a limited accuracy 
in the pure-state expansion or when working with (single- 
particle) scattering state o 83 ' 137 because these are not 
square renormalizable. However, the reservations apply 
equally well for MB collision problems (when there is no 
work around) treated in the formal LS approach] 28 ' 110 ' 138 
An evaluation of the Ehrenfest theorem becomes inde- 
terminate for open-boundary models which, by nature, 
must have infinite extension to avoid carrier depletion. 
Ehrenfest proof (for single-particle dynamics) rests on 
using Green's theorem 

/ (uAv - vAu)dr = / (uWv - vWu) ■ dr (8) 
Jv Js 

and on using an assumed cancellation of the surface 
terms (in 'S') at infinity^ 7 - There is no argument di- 
rectly available for a cancellation for scattering states 
(except, perhaps, by the single-particle work around sug- 
gested in Refs. 52). In fact, a careful handling of sur- 
face terms is essential for both the formal single-particle 
collision formulation,^ 7 - for a consistent NEGF evalua- 
tion of noninteracting tunneling transport ) 77 ' 139 and for 
a LS collision description of an open MB interaction 
problem j 109 ' 138 No argument exists for Ehrenfest theo- 
rem for the electron QKA which treats NEQ tunneling 
as a MB collision problem] 28 ' 62 ' 63 

Second, the inherent (GCE) thermodynamic charac- 
ter of tunneling 6 - 6 - and insight from GCE investigations 
of surface-physics interactions problems^^ 3 - (as well 
as from the exact steady-state ReNQSM^ 7 -^^^) give 
strong arguments for refining the NEQ description be- 
yond electrostatic forces. The insight identifies a need 
to improve the description of the Gibbs free energy ef- 
fects (as mandated by the Friedel sum rule) and to en- 
sure an adiabatic behavior in the implementation of nu- 
clei dynamics. Most tunneling systems can be expected 
to approach a steady-state behavio r 25 ' 95 ' 96 due to the 
open-boundary conditions and due to the presence of an 
entropy Row^^^l The QKA naturally treats such 
a system as an effective ensemble of tunneling systems 
and the QKA correctly reflects the existence of statistical 
fluctuations in the electron distribution at relevant times, 
i»iO. The fluctuations manifest themselves, for exam- 
ple, as current noise i 80 ' 140 ' 141 The presence of fluctua- 
tions al so m otivates a time-averaged operator approach, 
Ref . I93lll04l. 

Third, electrostatic forces will only be interesting if the 
assumption ([5]) provides an evaluation which is represen- 
tative for the ensemble. However, the source of an addi- 
tional force contribution with a statistical origin (for ex- 
ample, the electron-redistribution effects arising in open 
system o 67 ' 68 ' 92 ' 101 ) is also directly evident in the exact 
ReNQSM result for an effective (unnormalized) steady 
state NEQ density matrix* 6 - 7 - 

p hS = exp[-/3(H-Y hS )}. (9) 



This exact NEQ density matrix is given not only in terms 
of the Hamiltonian but also in terms of an operator-^ 2 - 
Yls f° r interaction-induced redistributions of electrons. 
Accordingly, it is important that the forces that effec- 
tively act on this ensemble (of steady-state NEQ tunnel- 
ing systems) reflect the sensitivity of this electron redis- 
tribution on the nuclei position J^ 3 - Moving to a consistent 
(exact) GCE thermodynamics theory is particularly rel- 
evant for combining an electron QKA with a NEQ BOA, 
for we seek a description of relaxation which mimics suit- 
ably defined adiabatic transformations. 

This paper takes a fundamentally different route than 
the Ehrenfest or electrostatic-force approach ^ to en- 
sure an improved description of an adiabatic character 
in relaxations and a better framework for implementing 
a NEQ BOA. Inspired by Hershfields exact reformula- 
tion of steady-state transport j 6 - 7 - this paper recasts a for- 
mally exact NEQ QSM theory, the QKA, in terms of a 
rigorous variational open-boundary NEQ thermodynam- 
ics theory. Being formally exact, it avoids approxima- 
tions which have previously been used to describe forces 
within a NEQ thermodynamical accoun t 51 ' 59 and I de- 
rive thermodynamic forces which are different from the 
electrostatic forces predicted and analysed in previous 
studie a 31 ' 48 ' 51 ' 52 ' 57 " — and which have gained wide usage 
in BTD implementations ! 40 ' 41 Because these forces arise 
from within a variational theory, the here-presented (ex- 
act) thermodynamical forces are given by a generalized 
Hellmann-Feynman (GHF) theorem, and they are explic- 
itly conservative. 



III. FORMAL THERMODYNAMICS RESULTS 

I use a continuum formulation to describe the NEQ 
electron dynamics in tunneling as a collision problem to 
ensure that the formal solution is nonempirical. That 
is, I emphasize a continuum description so that the for- 
mal NEQ thermodynamic solution permits calculations 
of parametric derivatives. The here-presented thermo- 
dynamic theory can (in connection with a NEQ BOA) 
therefore predict atomic-structure relaxations under op- 
erating (current carrying) conditions. Thus the present 
thermodynamic theory is, in principle, capable of deliver- 
ing an ab initio level of detail (as we can do in equilibrium 
with GS-DFT). 

This paper motivates and provides a definition of an 
exact NEQ thermodynamic grand potential value and a 
corresponding functional, 

ttcoiW = -ilnTr{e-^W-^W]}, (10) 
Q co i[p(t)} = Ttip^H^-Y^ + ^^K^W 

respectively. I denote the inverse temperature as /?. A 
subscript 'col' is used to emphasize that the NEQ in- 
teraction tunneling problem is solved as a formal MB 
collision problem for a general time dependence. This 
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approach was also taken in the formulation of the LSC- 
DFTi2& The functional formulation is a generalization of 
Mermin's equilibrium thermodynamic theory 1 ^ and it is 
closely related to ReNQSM descriptions. 67 ! 107 ! 122 ' 123 The 
NEQ thermodynamic grand potential (jlll) is a functional 
of the class of possible NEQ (time-dependent) density 
matrices p(t). The definition and formal results are made 
possible by the identification of an operator Y co \(t) for the 
NEQ Gibbs free energy for MB collision problems. 

This paper also demonstrates that the grand potential 
functional (|TT|) is extremal at the exact NEQ collision 
density matrix p co i(t)'- 



ficol [/>(*)] > ^col[j0col(i)] = ^col(i)- 



(12) 



That is, the grand potential functional acquires the ex- 
pected minimum (|10p as the extremal value and serves 
both to identify the exact NEQ grand potential value 
(HJJ]) and the exact NEQ density matrix for interacting 
tunneling problems. The Mermin-like variational prin- 
ciple ([12")) of this NEQ thermodynamic theory supple- 
ments the ReNQSM-based variational principle defined 
by a variance expression in Ref. [94| and supplements the 
GS, action, and MB-LS T- matrix variational principles 
used in DFT formulations. 

The variational property (|12p is proved in appendix A 
for members, p, of the class of normalized NEQ density 
matrices with a statistical weighting given by a Hcrmitian 
operator X(t) and normalization Ng/ t \, 



P e i N x(t)Pxlt)}' 
Px\t) = ex P[-/^(*)]> 



N 



X(t) 



= Ty {Pxlt)}- 



(13) 
(14) 
(15) 



The set of all possible NEQ solution density matrices, 
Pcoi(t) (resulting in MB collision problems) are contained 
in {-^x( t )/5x( t )}; the conservative nature of the QKA is, 
of course, sufficient to ensure that the normalization is 
independent of time. 

This paper uses the formal electron QKA of thermo- 
dynamic to identify new thermodynamic forces which I 
find explicitly satisfy the adiabatic criteria of a mening- 
ful NEQ BOA. Here I focus on tunneling systems in 
cases where a steady-state transport does emerge, as is 
expected to happen in general because of the entropy 
flow. 28 ' 95 ' 96 For these systems it is possible to define the 
solution grand potential as a state- function, f^LS- The 
thermodynamic forces are defined as minus the paramet- 
ric derivatives (with respect to nuclei positions Rj) of 
this state function, 



r GCE _ 

LS,R; = 



a 

9R 



LS . 



(16) 



These thermodynamic forces express a tendency of the 
(overall, electron-nuclei) system to maximize the (elec- 
tron) entropy and the electron redistribution (Gibbs free 
energy) by exploiting atomic-coordinate changes (thus 



extending entropy optimization beyond the maximiza- 
tion which already occurs at frozen coordinates) . 

I show that the thermodynamic forces are different 
from the electrostatic forces analyzed or used in a num- 
ber of previous studies i 40 i 48 i 51 ~ 54 i 57 ~ — The here-defined 
thermodynamic forces are given by a GHF variational 
theorem 



3 ft,-- a 



<9R 



on, 



H — Yls 



(17) 



The thermodynamic forces are explicitly conservative 
and I show that they are adiabatic. In summary, I find 
that use of thermodynamic forces is guided by a funda- 
mental principle and that they are the best approxima- 
tion for implementing a meaningful NEQ BOA. 

Finally, I note that the set of formal results leads to- 
wards a rigorous single-particle formulation of steady- 
state NEQ thermodynamics, a formulation that uses the 
LSC-DFT proof for uniqueness of density^ 



IV. PARTITION SCHEME AND 
HAMILTONIAN 

The partition scheme of Caroli et al, Ref. [H^, al- 
lows a description of tunneling which is first principle 
and predictive (i.e., parameter-free, atomistic). This is 
true when the Caroli partition scheme is extended from 
the original tight-binding framework to the continuum 
formulationi 28 i 63 The emphasis on an atomistic approach, 
i.e., with the description of the electron dynamics having 
a parametric dependence of the nuclei positions, is impor- 
tant because we desire a thermodynamic theory specific 
to the materials nature. Avoiding tight-binding parame- 
ters implies that electron behavior is completely specified 
by the position of the nuclei and by whatever external 
field that one chooses to apply. Subject to a suitable 
NEQ implementation of the BOA, one can therefore seek 
to relax atomic structure and morphology to provide a 
parameter-free and fully predictive thermodynamic de- 
scription of interacting NEQ tunneling (starting from an 
electron Hamiltonian). 



A. Initially disconnected equilibria system 

I assume, for simplicity, that the tunneling structure 
only comprises a left (right) lead r < zl(> zr) plus a 
central tunneling region 'C' in zl < r < zr. The con- 
tinuum (and therefore atomistic) Caroli partition scheme 
expresses the total kinetic energy K of the systems as a 
sum, K = Kl + Ku + Kc + SK, where Kl,r : c has sup- 
ports which are confined to the respective, separate com- 
ponents (left, right, and center) and where the connecting 
tunneling term SK is assumed to vanish identically in the 
far past, at t — > — oo and to exist exclusively at the inter- 
faces Zi = L/R between the partition-scheme components. 
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Using atomic units throughout, the total kinetic energy 
can be expressed by the operator 



The tunneling term in the continuum Caroli partition 
scheme can therefore also be written 



K = ' dv ' J *$(rO~<r'|fc>)V; s (rXl8) 6K = J2 J dr ^t( r ) 



-^ 2 + Mr) 



Mr) 



(r'|P|r) = -/in^r'-rO^-r). (19) 

In this second-quantization form one can easily define the 
tunneling term 

pZi+Sz/2 

SK = lim / dr x 

i=L/R a ->*t-oz/2 



I 



dr>t( r ')_ (r '|P| r) 4 (r) (20) 



The more traditional second-quantization form of the to- 
tal kinetic energy, 

K = Y,j (4 V2 ) Ur). (21) 

can also be used to obtain (below) an alternative formu- 
lation of the tunneling term SK. One is free to make ar- 
bitrary choices for the positions Zi—L,R, since the formal 
effects of the pertubation in 5 K (and all other interaction 
terms) are included to all orders. 

The continuum (and thus atomistic) Caroli partition 
scheme assumes that each of the disconnected subsections 
Hd = J2i=L c R-Hi is in equilibrium at different chemical 
potentials (J-l/c/r- Initially, at t — > — oo, the system is 
described by a static potential Vd(r) and corresponding 
operator 



V d = drv d (r)n(r). 



(22) 



Here n(r) = J2 S V'l( r )V ; s(r) denotes the electron-density 
operator. It may be convenient to assume that the form 
of t>d(r) for r reduces to a uniform background potential 
(j>L/R for z +/ — °o] m this case the values of <f>L and 4>r 
must be chosen to ensure charge neutrality of the initial, 
disconnected Hamiltonian at t — > — oo (and at z — > ±oo 
in general). The initial system in the continuum (and 
hence atomistic) Caroli scheme is simply 



H d = K i 

i=L/R/C 



V d . 



(23) 



I use {c\ i=L / R / c } to denote the set of operators which 
creates an electron in a H d SP eigenstate (of energy ex,i) 
with basis in the left, right, or center regions, respec- 
tively. The quadratic Hamiltonian (|23p can thus be ex- 
pressed 



E 



(24) 



i=L/R/C A 



(25) 



a form which perhaps is closer to the traditional start- 
ing point in NEGF calculations than (|20|) . The precise 
formulation of the tunneling term is of no consequence 
for the formal analysis presented in Ref. HH and here. 
The precise form is also irrelevant for the exact calcula- 
tions and for the effective single-particle computational 
approach which is also presented and developed here. 

A simple Gibbs-weighting operator X d and associated 
unnormalized density matrix, 



p d = exp[-ftX d ], 



(26) 



characterize the disconnected (t —¥ — oo) system at in- 
verse temperature ft. The initial density matrix is a di- 
rect product of L/C/R equilibrium density matrices and 
the initial Gibbs weighting can therefore directly be ex- 
pressed as a difference 



Xd — H d — Yd 



(27) 



between Hd and an operator for the t — > — oo Gibbs free 
energy, 



-L/R/C 



(28) 



The disconnected-equilibria (t — > — oo) density matrix 



exp[-ft(H d - Y d )} 



Tr{pd} Tr{exp[-ft(H d - Y d )}} 



(29) 



i=L/R/C A 



is independent of the value of the applied bias 0bias = 
<t>L — 4>R = — Hr- I note that the interacting tunnel- 
ing system will evolve to gain the correct density distri- 
bution and the relation 4>l — 4>r. = Mi ~ V>R i s neither 
generally valid nor used at relevant times t ~ 0. The 
initial density matrix depends exclusively on the initial 
electron concentration in the leads and the pair of values 
for lead chemical potentials, and, formally, on the initial 
occupation of the central island 'C'^ 



B. The interacting tunneling system 

The continuum Caroli partition scheme assumes an 
adiabatic turning on of the fixed tunneling SK, of a given 
MB interaction W (for example, the full electron-electron 
interaction) and of a static electron-scattering potential 
v sc (r) which includes the effects of the applied bias and of 
the set of atomic potentials. The scheme can also allow 
for an additional time-dependent potential </> g (r,i >t g ), 
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then describing a gate operation starting at some finite 
time t g . The total (time-dependent) collision potential is 



Ucoi(r,*) = fsc(r) + g (r,t). 



(30) 



This single-particle potential defines a quadratic contri- 
bution 



Vcoi(t) 



J dri> co i(r,i) n(r) (31) 



to the Hamiltonian. 

In several cases in the presentation of a thermodynamic 
theory for tunneling, the focus will be on the steady- 
state regime. Here the collision term simply describes 
the effect of the time-independent scattering potential 



/ 



drv sc (r)n(r). 



(32) 



In such steady state cases I will use a subscript 'sc' (or 
'LS') on formal evaluations to emphasize that the results 
are obtained when operators characterize the steady- 
state situation (or that the operators and expectation 
values are given by a formal LS solution). 

For a continuum description (which allows atomistic 
calculations) of electron dynamics in NEQ interacting 
tunneling, I investigate the full quantum-kinetics of the 
electrons in the connected (time-dependent and interact- 
ing) Hamiltonian 

H(t) = H d + a v (t) {6K + W + V col (t) - V d }, (33) 

The adiabatic turn on is described by a factor a v (t) which 
for all practical purpose*^ can be taken as exp(?7t), rj — > 
0; the original LS analysis 1 ^ of the MB collision prob- 
lems used exp(— r]\t\),T) = + . 

Furthermore, to aid the discussion of an exact, effec- 
tive SP formulation (for efficient calculations of the NEQ 
density, thermodynamics grand potential, and forces) in 
steady state, I also introduce the quadratic Hamiltonian 



=H d -V d + SK + V, 



(34) 



where the potential term V is described by a time- 
independent scattering potential, V = V sc . 



Ref. 11241 provides an excellent introduction and discus- 
sion of application of formal collision theory to condensed 
matter theory problems. The trick of a soft turn on of 
the perturbations was soon applied also for (infinite) con- 
densed matter problems to drive the system from a well- 
defined initial configuration and to the physically relevant 
case of an interacting MB system. The adiabatic turn on 
(and hence the inherent collision nature) is essential for 
the quantum-kinetics account of NEQ dynamics 1 ^ as well 
as specific NEGF calculations of interacting tunneling for 
it provides such calculations with a robust starting point. 

The identification of operators, Y co \(t) or Yls, for the 
emerging Gibbs free energy provides a link between the 
exact ReNQSM approach, formal MB collision theory 
(including LSC-DFT), and the present exact thermody- 
namics. The steady-state Gibbs free energy operator co- 
incides with the operator (originally termed l Y') that the 
steady-state ReNQSM uses to describe electron redistri- 
bution in the NEQ density matrix ©. This redistribu- 
tion operator can be interprete d 67 ' 68 



(35) 



X,i=L/R/C 



in terms of a family of operators, {^ LS A f } for MB- 
LS single-quasiparticle excitations of energy £x,i=L/R/C- 
These operators express the adiabatic evolution (de- 
scribed by a MB Liouville operator ChA = [H,A]) 
of the initial (independent) states subject to full MB 
interaction^iSSJM 



17] 



LS,X,i 



C 



H 



17] 



(36) 



This operator-based ReNQSM analysis can be general- 
ized and it is possible to express not only the general 
Gibbs free energy operator Y co \(t) but also the variational 
LS expression for the MB T-matrix^ 8 - in terms of single- 
quasi-particle excitations. A more complete account will 
be given in a forthcoming paper j 1 ^ below I simply sum- 
marize the findings relevant for the NEQ thermodynamic 
theory. 



NONEQUILIBRIUM QUANTUM 
STATISTICAL MECHANICS 



C. Relation to the quantum kinetics account 

I note that both the original Caroli scheme and the 
presently used continuum Caroli partition schem o 28 ' 63 
define a MB collision problem 109 for NEQ tunneling. 
Traditional MB collision problems involve an approach of 
originally free particles to a central interacting collision 
region; this approach can conveniently be expressed by an 
adiabatic turn on, for example, as used in the formal LS 
formulation? 1 ^ and in subsequent analysis by Pirenne^^ 
Dyson, 146 Gellmann and Goldbergerfii 7 - DeWittr^ 8 - and 
others. 



The formal structure of the exact NEQ thermodynamic 
theory is obtained in the spirit of Mermin's 1965 formula- 
tion of equilibrium thermodynamics^ In fact, in the ab- 
sence of any applied bias one regains the GCE system for 
which Mcrmin succeeded in formulating an equilibrium- 
thermodynamic DFT. 

A. Collision nature of transport 

The general strategy for studying MB collision prob- 
lems is to formally solve for projections of the full time- 
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evolution operator 

SW(t,-oo) 
1 dT~ 



H(i)U(t, -oo) 



(37) 



This operator describes, for example, the evolutions of a 
MB collision state 



|¥ f (t)) =W(t,-oo)|* e ) 



(38) 



originating and evolving from an initial state |$$). I shall 
assume that such initial states are simultaneous eigen- 
states of the disconnected-system operators Hd and Yd 
(of eigenvalues E^ and Y^, respectively) and I expand 
the unnormalized, initial density matrix 



Pd 



5>c> 



exp[-/3(^-r c )](^| 



(39) 



For a steady state collision problem one can explicitly 
construct a NEQ density matrix which incorporates the 
full complexity of the MB LS solution^ This LS solu- 
tion 



I* 



(+)\ 



IT] 



Ec — H 



IT] 



(40) 



constitutes an explicit construction of an interacting 
eigenstate of energy for steady-state problems. Steady 
state solutions arise in the collision problems when the 
collision potential reduces to the time-independent form 
(f3"2"j) and when the system dynamics is completely speci- 
fied by H(t = 0). The MB LS solution (gOj) constitutes a 
projection of the time-evolution operator U(t, — oo), since 
it satisfies 



exp(-iE^)|*[ +) ) 



= M(t,-oo)|* c > 



(41) 



in the limit of an infinitely slow adiabatic turn on, 
r\ -» + . Combining Eqs. and (|4"Tj) yields an (un- 
normalized) steady-state NEQ collision density matrix 



PhS 



M +) )cx P [-/3(^-y c )](* 



(+)l 



(42) 



The explicit construction (|4"2")l of the steady-state NEQ 
density matrix />ls provides a direct link of the NEQ ther- 
modynamic theory (and of ReNQSM) to formal MB LS 
solutions A closely related link was expressed already 
in the steady-state ReNQSM papers, Ref. l67l|68l as dis- 
cussed in subsection IV. C. 



B. Collision picture of interacting tunneling 

Standard QKA methods, including NEGF techniques, 
permit calculations of operator expectation values under 
NEQ conditions and must therefore provide a consistent 
treatment of the collision nature (|38p. The tested versa- 
tility of the MB LS equation for steady state problems 



motivates an effort to keep a corresponding collision pic- 
ture for solving general (time-dependent or steady-state) 
interacting tunneling through a NEQ collision density 
matrix p co \(t). 

Traditionally, the formal starting point for calculating 
the expectation value of an operator 0(t) in a NEQ sys- 



tem 1! 
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6 Ti{pdd n (t)} 



(43) 



The subscript 'col' emphasizes that a Caroli collision pic- 
ture underpins the definition. The evaluation (|43[) rests 
on the Heisenberg picture which represents a canonical 
transformation of operators, 

6(t) — > 6 H (t) = l/t(t, -oo)6(t)U(t, -oo). (44) 

The Heisenberg picture has advantages for time- 
independent problems (for it leaves the MB eigenstates 
fixed), but it is not a natural description in the pres- 
ence of the adiabatic turn found in the collision prob- 
lems. Nevertheless, the original t — > — oo system (Hd) 
has a well-defined separated-equilibrium thermodynam- 
ics pd and the formulation (|43[) still serves to uniquely de- 
fine the expectation value because the QKAs emphasize 
accurate (conserving) descriptions of the time-evaluation 
U(t,oo). 

For development of a NEQ thermodynamic theory it 
is more advantageous to make the collision nature (|3"8")l 
explicit in the operators for thermodynamic quantities. I 
work below with a collision picture defined by operator 
transformations 

i(f) — > A col (t) = U{t, -oo)A(t)U\t, -oo) (45) 

where A(t) denotes a time-dependent operator in the 
Schrodinger picture. This is done simply because it of- 
fers a simpler formulation and interpretation; the here- 
presented description is formally exact. This approach is 
suggested by the formal steady-state LS solutionis and 
it is here used for a general time-dependent tunneling 
case. In such a collision picture we obtain an explicit 
construction of the exact NEQ collision density matrix 

pcol(*) = U(t, -oo)p d U\t, -oo) 

= !**(*)> expH3(^-^)K*«(*)l- ( 46 ) 



The time evaluation of the density matrix is 



.rfPcol 

dt 



[H{t),pUt)\- 



(47) 



The formal time-evolution of other operators are also 
simpler than in the Heisenberg picture because of the 
opposite ordering of time-evolution operators in (|46[) . In 
using the collision picture for thermodynamic operators 
(1431) care is taken to verify that the resulting formulation 
remains exact and fully equivalent with the traditional 
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formulations of the QKA. 10 For the collision density ma- 
trix (|46p the formal equivalence 

/n\ ^- ^{PdOnjt)} Tr{p col (t)6(t)} 

{U/colW = rp r - i - 7^T~~ TTyj ■ I 48 ) 

Tr{p d \ Tr{p col (t)} 
follows simply from a cyclic permutation of operators. 



which obviously provides an explicit formulation of the 
Gibbs free energy in the original, disconnected-equilibria 
system. I note that, irrespectively of the interpretation, 
the operator Y co \(t) is uniquely specified by the collision 
nature of the Gibbs weighting operator X co \(t) or, equiv- 
alent^, of the exact collision density matrix (l5Tj) . 



Thermodynamical description of 
nonequilibrium tunneling 



The collision density matrix (|46j) can be expressed in 
terms of an evolving Gibbs- weighting factor X co \(t). I 
start from the disconnected-equilibrium Gibbs weighting 



H d -Y d 



(49) 



X co \(t -co) = X d 

and express the evolution 

X col = U(t,-oo)X d U\t,-oo), (50) 
0coi(f) = exp[-/3l col (t)], (51) 

thus identifying (|50l) as the emerging Gibbs weighting. 

The introduction of a NEQ Gibbs weighting factor mo- 
tivates, in turn, a definition of the thermodynamic grand 
potential 



^coi(i) 



lnTr{p co i(t)} 
lnTr{exp[-/3l col (<)]}. 



(52) 



The exact ReNQSM£2i£&22 shows that H - Y LS serves 
as a Gibbs weighting factor in the steady state problem. 
For the general tunneling problem I define 



Y col (t) = H(t) - Xcoi(t), 



(53) 



corresponding to a reformulation the normalized collision 
density matrix 



Pcolft) 

Tr{p col (t)} 



= expppJcoifi) - H(t) + Y co} {t)}. (54) 



A formal evaluation of the entropy S co i(t) in the exact 
collision density matrix can therefore be expressed^ 

«5coi(*) = -(ln(y3 co i)) co i + (ln(Tr{/5 co i})) co i 

= -fin^t) + P(H(t)) col - /3(Y col (t)) col (55) 

I consequently make the natural identification of Y co i(t) 
as the operator for the NEQ Gibbs free energy and I 
stress that the operator can be seen as emerging, with 
a time-evolution given by the LS collision picture (|45[) . 
This interpretation is consistent with the boundary con- 
dition, 

Y coi (t -> -co) =Y d = '^2 Pic\^c\ :i , (56) 



D. Nonequilibrium Gibbs free energy 

An actual evaluation of the exact collision density ma- 
trix (|51j) requires an explicit determination of the oper- 
ator for NEQ Gibbs free energy, Y co \(t). This was (with 
the interpretation stated above) a central achievement 
of the exact steady-state ReNQSM, proceeding through 
an explicit construction Yls order by order in the formal 
perturbation term^l Here I start instead with the colli- 
sion nature of the problem, to obtain a generalization to 
general tunneling and to clarify the explicit connection 
to a full NEQ thermodynamic theory. 

By the emerging nature of the Gibbs free energy oper- 
ator (|53l) . we can directly extract formal time evolution 



■ dY col (t) ~ dH(t) 

i — = [H(t),Y co i(t)\ + i- 



dt 



dt 



(57) 



a differential equation which, of course, must be solved 
subject to the (collision) boundary condition Y co \(t — > 
— oo ) = Y d . For the set of operators {c^ i=L / R / c } which 
(at t — > — oo) describes creation of single-particle H d 
eigenstates (of single-particle energies £\,i=L/R/Ci) we 
can formally establish the time-evolution in the collision 
picture 



M i (t)=U(t,-oo)c[ -Oo) 



(58) 



In turn, the formal solution ([55} provides an explicit de- 
termination of the operator for the emerging Gibbs free 
energy operator 



^col(t) 



J2 «$,i(*)lki(t). (59) 



X,i=L/R/C 



The exact NEQ Gibbs free energy value is calculated 
from (|59")) and from the uniquely specified NEQ solution 
density matrix, 



Y col (t) 



Tr{y co i(t)/3 co i(£)} 
Tr{p co i(t)} 



(60) 



For problems with a time-independent scattering po- 
tential v sc (r) , and when an actual steady-state transport 
results, the emerging Gibbs free energy Y sc must coincide 
exactly with the ReNQSM redistribution operator, 

^sc = Y LS . (61) 

The time evaluation is formally described by 

Y sc (t) = U(t, -oo)Y d U\t, -oo), (62) 
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but it is alone the MB Hamiltonian [H = H(t = 0)] 
and the open boundary conditions which determine the 
properties of this emerging Gibbs free energy. The formal 
time dependence is therefore limited to the adiabatic turn 
on itself, 

Y sc ~Y d + exp(r)t)[Y BC (t = 0) - Y d ], (63) 

and there must be a cancellation of the time-dependence 
in the creation and destruction operators in (1621) . As a 
consequence, the resulting Gibbs free energy is charac- 
terized by a time evolution 

which is identical to that derived for the steady-state 
electron-redistribution operator Yls derived and analy- 
sized by Hershfield and S chiller .^I^ 



energy content. The appearent difference arises be- 
cause the formal MB LS solution to the collision prob- 
lem implies a projection (further discussed in appendix 
B) of the full dynamics upon the eigenstates of the 
(i w 0) Hamiltonian. 109 ' 124 ' 147 However, for extended 
systems it is then essential to also include effects of 
both collision-state renormalizatio n 124 ' 148 (which arises 
whenever the problem is not asymptotically free) and of 
level shift o 145 ' 147 (which, for example, reflects the effects 
of band-structure formations on quasi-particle energies 
e x ,i). Since the steady-state ReNEQSM, the LSC-DFT, 
and the here-presented thermodynamic theory make use 
of the MB LS projection for steady state tunneling, there 
is no relation between the normalization of the H d and 
the H eigenstates, and there is no relation between the 
normalizations of the disconnected-equilibria and NEQ 
solution density matrices. For condensed matter systems, 
we have no binding of the steady-state total internal en- 
ergy content (|66[) to an arbitrary initial energy spectra. 



E. Total-internal energy and state renormalization 

When dealing with a fully time-dependent problem we 
obtain the exact total internal energy content directly 
from the NEQ density-matrix solution 



U co i(t) 



Tr{Hp col {t)} 
Tr{p co i(t)} 



(65) 



This expression has just a trivial denominator, 
Tr{/5 co i(t)} = Tr{pd}, because of the unitary character 
of the time-evolution operator, and the conserving na- 
ture of the NEQ solution density matrix p co i(t)- 

At the same time, it is important to emphasize an ap- 
pearent (but not actual) difference which arises in the 
evaluation of steady-state thermodynamic quantities de- 
pending on whether one retains the full NEQ time depen- 
dence or whether one exploits the ReNEQSM [and the 
connection to formal MB LS solution via (|42p]. The ap- 
pearent differences are most clearly evident in the dcfin- 
tion of a meaningful, i.e., system-specific evaluation of 
the steady-state total internal-energy content, 



Tr{/3 L s} 



(66) 



Because the solution density matrix /Sls commutes^ with 
the fully interacting Hamiltonian (H ) , it is tempting (but 
wrong) to conclude that the set of initial t —> -co dis- 
connected equilibria eigenstates (along with the initial 
thermal distribution) provides a complete specification of 
C/ls- Such a conclusion would not only render the ther- 
modynamic description of total internal energies mean- 
ingless but would also be in conflict with the formally 
exact, system-specific evaluation (|65l) . 

For a steady-state tunneling problem there can be 
no difference between the formally exact evaluation and 
the steady-state description. The evaluation must 
be a meaningful determination of the total internal 



VI. VARIATIONAL NONEQUILIBRIUM 
THERMODYNAMIC GRAND POTENTIAL 

Section III provides a summary of the central results of 
this paper, namely a description of the nature and prop- 
erties of the exact NEQ thermodynamic grand potential. 
Appendix A contains the proof of the variational prop- 
erties, generalizing the equilibrium results of Mermin. 16 
Here I only summarize the findings. 



A. Notation and motivation 

The evaluation of expectation values [for a general op- 
erator A(t)] is viewed as a functional 



(Mt))n 



Tr{p(t)A} 



(67) 



Tr{p(t)} 

of a general NEQ, time-dependent density matrix p(t). 

For a general NEQ density matrix (electron distribu- 
tion) there exists a natural functional evaluation of the 
entropy 



sim] 



(68) 



-(\np(t)) m + (H^{p(t)}]) m - 

Similarly, the evaluations of the NEQ internal energies 
and Gibbs free energy 

U col [p(t)} = (H(t)) m , (69) 

Y col [p(t)} = (Y co i(t)) m , (70) 

are set directly by the Hamiltonian H(t) and by the 
emerging Gibbs free energy operator Y co \{t). 

The NEQ evaluations of entropy, internal energy, and 
Gibb free energy functionals suggest, in turn, a functional 
form for the NEQ grand potential: 

ncci(t)[p(*)] - ~s[p(t)} + UcoiW)} - Y col [p(t)}- (7i) 
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This expression is a generalization (fTTj) of Mermin's 
(equilibrium) result to a NEQ thermodynamic grand po- 
tential for collision problems. The formulation ([7T|) is 
relevant for tunneling systems described by the contin- 
uum Caroli partition scheme. 

For a compact presentation of derivatives and forces 
it is convenient to view the NEQ thermodynamic grand 
potential functional (fTTj) as a special (collision) instance 
of the more general functional 

n A(t) b(t)] = (x(t) + /r 1 inp(t) - r 1 HTr{p(t)}}) m - 

(72) 

Obviously, we have 

^coi[p(i)] =%(t)=H(r)_r col (t)[p(*)]- ( 73 ) 



B. Variational thermodynamics 

The variational and extremal property, Eqs. (|12p. ap- 
plies among the class of density matrices which is defined 
in (|14[) . This class of density matrices obviously encom- 
passes all possible NEQ density matrices that emerge 
(with conserving dynamics) in a general collision prob- 
lem, Eqs. (|5D|) and (|5T|) . The extremal property pro- 
vides us with a tool to identify the correct NEQ collision 
density matrix f|46[) for a given (interacting) tunneling 
problem specified as continuum Caroli partition scheme, 
Sec. IV. 

As in Mermin's analysis one immediately finds a sta- 
tionary property of 



X(t) = ^x(t)[Px(t)\ 



(74) 



which can be formulated from the functional expressions 



an 



X (i)+AA(t) 



dX 



A=0 



\ Q X Px (t)+XA(t) 



A=0 



(75) 



As in the corresponding equilibrium case,— this property 
is central to demonstrating the variational and extremal 
properties of the thermodynamic grand potential at the 
solution density matrix p co \(t). 



C. Extremal property as a tempered maximization 
of entropy 

The NEQ solution density matrix p co \(t) is unique and 
can be identified by the minimum of the variational thcr- 
modynamical grand potential functional. From the con- 
struction of the thermodynamic potential it follows that 



15_ 

Wp 



s[p] 



Sp 



{n col (t)[p} - (u col (t)[ P ] - Y col (t)[p})} , 

(76) 



and that, following a system perturbation, p{t) ^ /3 co i(i), 
we expect a relaxation because it corresponds to a max- 
imization of electron entropy. This maximization is 
naturally limited or tempered by a functional for the 
total-internal-energy and Gibbs-free energy difference, 
U col (t)[p]-Y col (t)[p). 

The maximum-entropy character of the here-presented 
variational thermodynamic theory makes it clear that 
there are very close relations to the range of constrained 
maximum entropy formulationsi 73 ' 74 i 76 However, the are 
also important differences. While the contrained max- 
imum entropy approaches use Langrange multipliers to 
specify selected expectation values (for example, the cur- 
rent flow), the explicit need for constraints is bypassed 
in the ReNQSM and in this variational thermodynamic 
description. Specifically, I stress that the approach does 
not impose a value for the internal energy or Gibbs free 
energy content and that U co i(t)[p] — Y co i(t)[p] is itself 
uniquely specified by the boundary conditions and by 
the actual time evolution. 



D. Generalized Hellmann-Feynman theorem 

The fact that the exact NEQ solution can be viewed 
as a maximization of entropy (subject to exact evalua- 
tions of thermodynamic boundary conditions) also pro- 
vides a rationale for introducing generalized thermody- 
namic forces 

d_ 

dp" 

The more specific thcrmodynamic-forcc definition 



F c G ol^N-^c ol (t)[p] 



(77) 



d 



F^(*) = -0jrfi=oi(*) 



(78) 



is given by a parametric dependence of the density matrix 
(as will be further motivated for steady-state problems in 
the following section) . The combined electron- nuclei sys- 
tem is here simply assumed to further the maximization 
of the entropy content (in the electron system) . 

A GHF theorem exists to simplify the determination 
of such thermodynamic forces, here as discussed for time- 
dependent tunneling. As an explicit expression of (|75p . 
one can specify the perturbation 

A W = «5-^coi(t) ■ SILi (79) 

and evaluate the changes corresponding to a derivative 
of the thermodynamic potential along a coordinate-line 
element (5EL: 



d 



dX c 



This provides a demonstration of a NEQ GHF, 



i) 



^coi(i) 



a 



H - Y coX {t) 



'Pcol(t)' 



for a thermodynamic description of tunneling. 



(80) 



(81) 
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E. Noninteracting transport, link to the 
single-particle Lippmann-Schwinger solution 

For noninteracting particles, described by the 



quadratic Hamiltonian Hy , 
SP LS solutions^ 



it is clear that the set of 



LS,A,i/ 



IT] 



C H (0) 



11] 



(82) 
(83) 



constitutes the basis of what is formally the MB solu- 
tion. These states and creation operators retain a gen- 
uine independent-particle nature and hence the emerging 
Gibbs free energy, denoted Yy is also quadratic. It is 
also clear that these SP solutions must enter both the 
NEQ solution density matrix and the thermodynamics 
value. In particular, it is possible to obtain a separable 
form of the NEQ solution density matrix 



TT «(°) 
"V-/'; s .\.. 

exp[-/3(e A ,i - Hi)h\. 



(84) 
(85) 
(86) 



and of the corresponding NEQ thermodynamical grand 
potential. It is instructive to present a formal argu- 
ment that one can derive the set of formal SP LS so- 
lutions, as a direct consequence of the here-presented 
variational thermodynamical description (for a quadratic 
time- independent Hamiltonian). 

It is sufficient to limit the variational space (for the 
NEQ solution density matrix) to that produced by a 
corresponding set of mutually orthogonal trial functions 
V'A.^r) = (r|^A,i), each of these approximating the form 
a SP scattering state from Hd eigenstates of energy ex,i- 
No assumption is made for the normalization of these 

trial scattering wavefunctions whereas I use 4 } to de- 
note creation operators for such states with a suitable 
standard (delta-function) normalization 



\i>> 



(^X,i\^. 



X.i 



4/2 



(87) 



I also introduce n\j — (i>\i and define the trial den- 
sity matrix solution 



P 



Pa,-. 



H-X,iPX,i, 

exp[-/3(e A , 4 - (M)n\,i\, 

= l + hx ti (exp[-P{ex,i-Hi)]-l). (89) 

I note that the difference between the internal and Gibbs 
free energies is a SP operator which can be expressed in 
terms of the trial wavefunctions 



(o) o(o) 



EE 

A' ,i* X.i 



(Oh 



(90) 



The variation form of the grand potential value is ob- 
viously also separable 



(0) 
V,X,i 



(91) 
(92) 



(o) 

V,X,i 



ilnTr{^ A;i } 



--ln(l + exp[-/3(e A , J - Ml )]) (93) 



Tr{p A;i [^ 0) - + (1//?) HPx^ ± ) 



K0) 



Tr {PA,i} 



I observe that the components of the grand-canonical po- 
tential functional has the minimum 



^v,X,i\Px,i = phS,X,i] = ®>v,\,i- 



(95) 



Below, in the spirit of Ref. Il49l I investigate the SP ther- 
modynamics by considering the variational contributions 
from each of the orthogonal trial scattering eigenstates, 
{A, i}, at a time. 

It is straightforward to evaluate the terms in the vari- 
ational separable grand potential functional, 

.ilPl = («,< - W)(*M>fc., (Cfe J - l ) ( 96 > 



P\,i 



(ex,i - l^i)(ipx,i\i>x,i) 
Tr{px,jnx,i} 
Tr{pxA 

exp[-^(e A , t - Hi)] 
1 + exp[-/3(ex,i - Hi)] 



(98) 



It follows that the variational nature of the NEQ ther- 
modynamical grand potential can be expressed 



^A,iW 



^<LL5a, 
^A,i(r) 







(99) 



;>[ i-l^V,A,*[^A,i] i/rr(0) t>(°)m 7 \ 

D^xA = ( r \i H v ~ Y v Mx,i) 



Dh[>x 



C[i>x,i](ex,i - Hi)(T$titi0) 
(101) 



Since C[ipx,i] — 1 (at the solution), it is possible to ex- 
press the variational condition (|100l) simply as the con- 
dition 



(ijf -Yy">)tyx,i) = (ex,i-HiMx,i). 



(102) 



This equation is obviously solved by the formal SP LS 
solution, (1521 and 
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F. Relation to standard ballistic transport codes 

It is also interesting to discuss the nature of the widely 
used BTD approach for calculations of the electron den- 
sity (and transmission properties) in NEQ tunneling i£2r— 

Again, I consider only steady-state transport in a 
time-independent scattering potential and introduce 
both a system-specific potential-energy and (mean-field 
Coulomb) Hartree interaction terms 



V[n] = J 



drv sc (r)n(r). 



= -Idr 



dr 



,n(r)n(r') 



(103) 



(104) 



Both are here expressed formally as a density functional, 
an approach for which I present formal arguments in the 
following section. Of central importance for the useful- 
ness of the BTD method is undoubtably also the empha- 
sis on using a density functional expression to approxi- 
mate effects of the MB interaction. The BTD chooses 
to approximate these by using the exchange-correlation 
functional term of GS-DFT, that is, approximations to 
the universal functional for equilibrium systems: 



E xc [n] 
G GS [n] 



G GS [n]-G<®[n] 
{9 \K + W\9 )-*s[n] 



(105) 
(106) 
(107) 



specified as a difference between expectation values of 
the fully interacting GS |\Po) and of the noninteracting 
GS |\&o°^). From these functional the BTD constructs an 
effective scattering potential 

UbtdM(i-) = v BC (r) 



5$ 



H 



Sn(r) Sn(r) 



(108) 



for ballistic transport, that is, SP scattering 
dynamics 

The fact that solutions of the BTD dynamics are ex- 
pressed in terms of SP LS solutions, ensures that the 
BTD method (for electron-density calculations) is in- 
herently consistent, automatically satifying conservation 
laws. The BTD method incorporates the essential Gibbs- 
free energy effects as it is consistent with the charge ad- 
justments stipulated by using Freidel sum rule on dynam- 
ics in the effective potential (|108p . In effect, the BTD 
method thus consists of constructing not only the effec- 
tive potential ubtd [n] but also of constructing and using 
the associated SP Gibbs free energy operator, 



/M'AVbtd.V'' WbtdA'J' 



(109) 



X.i 



from SP creation operators (ipy g TD , for the BTD SP 
LS solution, 

A determination of the Gibbs free energy operator and 
the Gibbs free energy effects on consistent GCE thermo- 
dynamics forces (as described at the BTD level of ap- 
proximation), is essentially already available in standard 



BTD codes since these codes already evaluate the scat- 
tering states (in the iterative BTD determination of the 
electron density). 



VII. DENSITY FUNCTIONAL THEORY OF 
INTERACTING TUNNELING 

The proof of the variational property (TT21 of the 
NEQ grand potential does not require, nor by itself pro- 
vide uniqueness of density in NEQ tunneling, Ref. [28]. 
This is in contrast to the case of thermal properties in 
equilibrium*^ The emerging (rather than fixed) expres- 
sion for the NEQ Gibbs free energy operator Y co i(t) does 
not seem to provide any reductio ad absurdum argument; 
changes in the term V co \ will make themselves felt both in 
the Hamiltonian and in the Gibbs free energy operator. 

Nevertheless, the previously formulated LSC-DFT still 
provides a proof for uniqueness of density 28 for collision 
problemsi^S and defines a SP formulation for formally 
exact calculations based on universal functionals. The 
formally exact theory contains the set of BTD descrip- 
tions as a first consistent but also unquantified approxi- 
mation. The here-presented thermodynamics theory pro- 
vides a framework for systematic improvements beyond 
BTD. 



A. Uniqueness of density in nonequilibrium 
thermodynamics 

The closely related LSC-DFT— provides an indepen- 
dent proof for the uniqueness of density and hence ob- 
tains a variational density functional account. The uni- 
tary evolution provides an explicit construction for the 
mapping 



A/" : v col (r,t) — ► n(r,t) = (n(r)) col (f). 



(110) 



The LSC-DFT proof for uniqueness of density, Ref. [28| ex- 
plores the time-evolution in the collision picture to prove 
that Af is unique, i.e., that for every v co i(r,t) there can 
be at most one time-dependent density. Also, the case of 
steady-state tunneling resulting from the adiabatic turn 
on is just a particular form of time-dependent behav- 
ior and the uniqueness of density extends also to cases 
when the collision term reduces to a time-independent 
scattering potential t> sc (r) . Under the assumption that a 
steady-state behavior emerges from v sc (r) the mapping 
Af is then instead completely described by the steady- 
state density matrix ©: 



Af : v sc {r) — ► n(r) = (n(r)) 6 



(111) 



It follows cither by the boundary conditions discussed in 
Ref. [67| or directly by differentiating the mapping (|110j) 
that the resulting density can have no more than an adi- 
abatic time dependence, n(r,t) ~ n(r) ex.p(j]t). 
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This unique mapping between density and external po- 
tential allows an explicit construction of density func- 
tionals for the thermodynamics of NEQ tunneling. I 
observe that it is possible to choose a detailed form of 
the adiabatic turn on such that it is compatible, even in 
a stringent mathematical sense, with both the require- 
ment of LSC-DFT2 8 . and with the present NEQ ther- 
modynamic description™; this is discussed further in 
Appendix B. A robust, universally applicable, linking of 
the LSC-DFT with the NEQ thermodynamics also re- 
quires an assumption of a V-representabilityj 1 ^ i.e., that 
for every NEQ electron density there exists an external 
potential for which the solution density matrix produces 
that density^ 

A solution NEQ density matrix corresponding to a 
given density n'(r, t) can therefore be viewed as a general 
density functional 



Pn'(t) = p[n'(r,t)] = Pcol,V'<-n'(*)- 



(112) 



Below I focus the discussion on steady-state transport in 
a time- independent scattering potential (|103[) . 

The formally exact description of the overall dynam- 
ics in NEQ tunneling also requires the introduction of a 
number of other density functionals 



F[n] 
Y v [n] 



Tt{p n [K + W + (l/j3) ln(p n /Tr{p n })]i 



Tr{/5„} 



Tr{/5„} 



113) 



(114) 



The variational character of the thermodynamical poten- 
tial can be expressed 



n v [n] = V[n] + F[n]-Y v [n] 

dn v \n] on 



(0)r 



On 



dn 



d_ 

dn 



OyFnj 



0, 



(115) 



(116) 



as a correction to the SP dynamics derived in Subsec- 
tion VI. E. As the subscripts l V indicate, however, the 
evaluation of Gibbs free energy term (and this formu- 
lation of the density functional for the thermodynamic 
grand potential) must involve an explicit construction of 
a (nonuniversal) operator Yy . 



shows that it (partly) contains the BTD method 8 - - — as 
a lowest-order approximation which stills retains a con- 
sistent account of Gibbs free-energy effects in the eval- 
uation of the electron density (but obviously not when 
BTD is subsequently used to determine relaxations from 
electrostatic force o 40 ' 85 ). 

I observe that since the scattering potential is a func- 
tional of the density, Eq. (|103l) . it follows that the LS 
operator solution can also be expressed as a density func- 
tional 



LS,X,il 



in 



C 



in 



(117) 



From this functional form of the MB LS operator solu- 
tion, it is possible density derivative of the variational 
thermodynamics grand potential in terms of universal 
functionals. I define 



Y[n] 
Y'[n] 



X^ i ^LS,A,i[ n ]V'LS,A,i[«] 
i 

Tr {p n Y n } 

Tr{/5„} 
Tv{p n {dY n /dn)} 
Tr{/5„} 



and note that 



£ Yv [n] = ±Y[n]-Y'[ n ] 



(118) 

(119) 
(120) 

(121) 



express the derivative of the Gibbs free energy in terms 
of universal functionals. 

We seek a SP formulation where solving for fictitious 
SP scattering wavefunctions ^> e ff,A,i( r )j moving in some 
effective scattering potential V c s, constitutes a complete 
specification of both the electron density and the thermo- 
dynamical potential fly (of the fully interacting system) . 
The individual SP scattering state is a SP LS solution 
which emerges from an Hd eigenstate A, i (of energy 
Assuming a suitable normalization of such solutions leads 
to a simple determination of the elecron density variation 



n(r)=^|V>eff,v( r )l' 



(Tlx 



l/Pn 



(122) 



X,i 



B. Universal functionals and formally-exact 
single-particle descriptions 

While the grand thermodynamical potential functional 
(|115[) contains a nonuniversal component (beyond V[n]), 
it is still possible to achieve an effective single-particle for- 
mulation in terms of universal density functionals. This is 
an important point for the variational thermodynamical 
theory does therefore enable an efficient evaluation for 
general interacting tunneling problems. Moreover, the 
following development of a formally exact SP formulation 



where {n\.i) p„ is evaluated in Sec. VI. E. To express a 
rigorous SP formulation I furthermore introduce a set of 
exchange-correlation terms 

F KC [n] = F[n] - F (0) [n} (123) 
^xcH = Y[n]-Y (0) [n] (124) 
Y x ' c [n] = Y'[n]-Y''W[n] (125) 

where, as before, superscripts '(0)' identify functional val- 
ues evaluated at W = 0. Using the framework from the 
formal derivation of the SP scattering states from the 
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variational potential, Sec. VI. E, it follows that the ficti- 
tious dynamics is uniquely specified by an effective equa- 
tion of motion 



(0) 



^ * ff ,A,(r) 



dn 



^eff,A,i( r ) 

SY m [n] 



(126) 



6n(r) Sn(r) 

Since we have assumed a suitable normalization, the 
equation of motion for an occupied SP LS solution wave- 
function can be reformulated, 



(r\(H { v 0) -Y^M cS ^) + 



(oh 



thermodynamical DFT. It is clear that the BTD method 
is consistent: it is built from SP LS scattering states 
and it follows that the description (in itself) is compati- 
ble with the Friedel sum rule. The BTD method captures 
Gibbs free energy effects [to the extent that these arise in 
the BTD choice of the effective potential (|108[) ] through 
the implicit construction of the quadratic operator Vy BTD , 
+ y xcN ) ^cffXi^) Eq. pUSjl . Within the BTD calculation method, where 
there exists a relation between the density and the orig- 
inal potential V, one can also define the BTD solution 
density matrix as a density functional /5^ TD and, in turn 
a density functional for the BTD approximation for the 
Gibbs free energy in actual potential V: 



Sn{r) xc[ J 5n(r) 

C[lpeff,\,i]{e\,i ~ M0( r l^eff,A, 



(rIV'eff,. 



(128) 



The real- valued functional V x ' c [n] can be viewed either as 
a constant potential shift or as a shift in the eigenvalue 
e\,i, for example, as expected when an infinite scattering 
potential causes a change in band-structurei 145 ' 147 

The SP LS scattering-wave solutions (for dynamics in 
an effective SP potential V e s that gives rise to the SP 
Gibbs- free energy operator Yy dt ) is: 



(r|(4^-yg)l^ 



(e A , l -M l )(r|Vcff,A, l ). (129) 



An assumption of V-representability is, of course, re- 
quired to conclude that such an effective potential V e s 
actually exists for a given solution electron density n(r) 
of the interacting problem. However, given such an as- 
sumption, it is clear that the terms of the fictitious SP 
fomulation (|129[) and the variational evalulation (|128[) 
must coincide and this provides a framework to identify 
the uniquely defined effective SP potential 



V eS [n] 



dr v eS [n}(r)n(r) 



(130) 



At the same time, and unlike in GS-DFT and in BTD, 
it is not so that a comparison of the dynamics specified 
by the variational principle leads directly to an explicit 
determination of a new effective potential, 



V eS [n] ^ V[n] + F x , 



(131) 



To search for an effective SP potential for description 
of MB interactions effects in NEQ tunneling we must 
instead require an implicit approach with concerted ad- 
justments of both the potential form and of the associ- 
ated Gibbs free energy term. Future works will explore 
approaches for simplifying this implicit search, as well 
as possible approaches to determine the set of universal 
functionals which are required for the here-suggested SP 
calculational scheme. 

Finally, we are now able to discuss the nature of BTD 
method for calculations of electron density as a lowest- 
order approximation to the above-discussed, exact GCE 



Y 



BTD r 



Tr{/3g TD ^ BTD } 
Tr{pBT D } ' 



(132) 



At the same time, it is clear that the BTD method does 
constitute an approximation in that it does not express 
the thermodynamics in terms of the actual MB function- 
als F xc [n], Fa; C [n], and 3^ c [n]. Comparing the expressions 
it is clear that a systematic improvement of the BTD 
method can be made by also consistently including the 
MB NEQ scattering effects which are expressed in the 
functional-derivative difference 

(f^m-Bb-M-^M)- 



Sn(r) 



5n(r) 



(y xc [n]-y, 



y BTD [n]+4%l)-+ 



F x >] ? 0.(133) 



VIII. STEADY-STATE THERMODYNAMICS 
AND ADIABATIC FORCES 

The formulation in terms of a NEQ grand potential 
functional (above) is by construction fully equivalent to 
the exact solution in the QKAi^ - Specifically, a mini- 
mization of the general NEQ grand potential functional 
difference f2 co i[p(t)] — tt co i(t) will uniquely identify the 
exact NEQ density matrix p co \(t). Access to a formally 
exact thermodynamic theory provides a precise starting 
point for a discussion of adiabatic relaxations by suitably 
defined thermodynamic forces in the NEQ tunneling sys- 
tems. 

Needless to say, a full discussion and evaluation of 
wind forces and of electromigration is beyond the present 
analysis. The same applies for non-adiabatic forces (fol- 
lowing sudden system changes) and for tunneling sys- 
tems where the electrons exhibit an unstable, intermit- 
tent behavior^ To emphasize that I am working with 
steady-state tunneling where the MB LS solutions en- 
sures an automatic satisfaction of all conservation laws, 
I choose here and in the rest of the paper to replace the 
subscripts 'col' by 'LS' on the thermodynamic expres- 
sions. While I could have continued to use the subscript 
£ V (as was done in the discussions of SP formulations 
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of steady-state results) , I am thus also emphasizing that 
the following discussion of thermodynamics forces is a 
result which is independent of the density-functional de- 
scription, Sec. VII. On the other hand, the BTD method 
or even a full density functional SP formulation will pro- 
vide an efficient evaluation of the here-presented force 
description. In any case, to simplify notation, I also as- 
sume below that the time-independent scattering poten- 
tial v sc (r) is entirely specified by the position of the fixed 
nuclei {R^}. 



A. State functions for steady-state tunneling 

I explore the NEQ thermodynamics of tunneling sub- 
ject to the assumption that the NEQ density matrix 
(and thus all physical properties) has in fact relaxed to 
a steady state. To create a framework for a discussion 
of a NEQ implementation of the BOA, I investigate the 
variation in the set of thermodynamic quantities with in- 
finitely slow deformations of the morphology (that is, in- 
finitely slow changes in the set of nuclei positions {Ri}) 
and I inquire what forces would then guide the relax- 
ations arising with the tunneling transport. 

An exact evaluation of the steady-state NEQ den- 
sity matrix pls directly provides the corresponding exact 
NEQ evaluations 

Uls(vl,VR, McS {Ri}) = U h s[phs], (134) 

Yls(hl, Hr, He; {Ri}) = Y hS [phs], (135) 

SLs(A*i)A*R 5 Mc;{Ri}) = 5[/3ls] (136) 

^LsCMij/^McsiRi}) = ^ls[/5ls] (137) 

for the internal energy, the Gibbs free energy, the entropy, 
and of the NEQ grand potential, respectively. The func- 
tions (|134l) - (|137[) are by construction state functions of 
the indicated parameters. 



B. Surface-state mediated adatom interactions on 
Cu(lll) 

As an illustration and to simplify a discussion (below) 
of the nature of forces in GCE thermodynamics, it is 
instructive to summarize details 1 ^ of a GCE thermo- 
dynamics stud y 101 ! 106 for an equilibrium surface-physics 
problemji2£ The GCE study of the mutual adatom 
interaction mediated by the metallic Cu(lll) surface 
state formally proceeds within Mermin's thermodynam- 
ics DFT by combining the Harris-scheme 1 ^ and a formal- 
scattering theoryi^ evaluation of the single-particle con- 
tribution to the interaction energy. 

The Cu(lll) adatom-adatom interaction problem is 
at the same time relatively simple and demanding. It 
is relatively simple because the isotropic surface state 
has an approximative free-electron dispersion 1 ^ given 
by effective mass m e s. It is also demanding since the 
surface-state scattering off the adatoms (and hence the 



mutual adatom interaction) must be described nonper- 
turbatively. Experimental observation of the Fermi-level 
local density of state oscillations 1 ^ shows that while the 
surface-state scattering off an adatom can be approxi- 
mated as having an s-wave character, the s-wave scat- 
tering phase-shift is very large, 5 F ~ ±7r/2. This non- 
pertubative surface-state scattering reflects the strong 
scattering-induced leakage of surface-state electrons into 
bulk states*!^ At the same time, the presence of the 
bulk-state electrons also provides a key simplification 
in the Harris-scheme evaluation i 1 1 1 1 5 We can use a 
general surface-physics expectatio n 101 ' 154 that charge- 
oscillations in the surface state is almost exactly cancelled 
by opposite oscillations by bulk-state electrons. It follows 
that it is sufficient to focus on evaluation on the single- 
particle terms in the GCE formal scattering theory* 1 ^ 

In the presentation of details in a GCE Harris-scheme 
calculation of mutual adatom interactions on Cu(lll), I 
assume that the bottom of the surface-state band is lo- 
cated at 0, while a Fermi energy e F identifies the chemical 
pote ntial in the zero-temperature evaluation. Following 
Ref. llOll I introduce Green functions 



TOeff 



g° R (d] e) = '-^-K (-iq e r + irj) 



(138) 



for the two-dimensional surface-state dynamics. Here 
e = h 2 q^/2m c s and Kq denotes the two-dimensional 
zeroth-order Bessel function. Since the experimental and 
theoretical characterizations rests on a s-wave description 
of the adatom-surface-state interaction, I also introduce 
the surface-state T-matrix 



tis(e) 



it? 



J2S (e) 



(139) 



for multiple (5-function scattering! 101 ' 104 This energy- 
dependent T-matrix is given by the s-wave phase shift 



So(e) — arccot [ln(e/ep)/7r + cot(<5p)] , 



(140) 



which, in turn, is entirely specified the experimentally 
observed Fermi-level (s-wave) phase shift, S F « ±7r/2. 

The Harris-scheme evaluation for the GCE formal scat- 
tering theory analysis focuses on the evaluation of the 
changes in density of states £>2s(e, d; S F ) which arise when 
two adatoms approach from an infinite distance and 
to some finite separation d. Using formal scattering 
theory 153 it is natural to express this change in terms 
of the integrated density-of-state variation 



V 2s (e,d;S F ) = ^-AN 2s (e,d;6 F ) 
oe 



(141) 



AN 2s (e,d;6 F ) = --Sin {l - [t ls (e)g R (d; e)} 2 }U2) 

In turn this allows a precise evaluation of single-particle 
contribution in the GCE thermodynamics description. 
In fact the Harris scheme permits the approximative 
1/(8 — !• (GCE) determination of the thermodynamics 
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properties 



by the variational property 



U&(d;6 F ) 

YH(d;S F ) 
nl s s (d:5 F ) 



F de'e'-^-jAN 2s (e\d;6 F ) (143) 



e F AV 2s (e F ,d;c5 F ) 
U? s (d;5 F )-Y£§(d;6 F ) 



(144) 

f) 

[ F de'AiV 2s (e',d;5 F ). (145) 
Jo 



The Harris-scheme and nonperturbativc formal- 
scattering-theory approximation for the thermody- 
namical grand potential, Eq. (|145|) was evaluated in 
Ref. llOll and soon after used for the comparison with 
surface-physics experiments in Ref. Il03l 



C. Electrostatic atomic forces 

It is possible to ignore all thermodynamics, that is, all 
effects of including the Gibbs free energy, when calculat- 
ing the forces on the nuclei 



-ics,NEQ 
LS,R; 



,dH_ 
v dR/ pLS ' 



(146) 



As above, the forces (|146|) are identified by a superscript 
'es' since they are purely Coulombic or electrostatic in 
nature. A number of paper o 48 i 52 ~ 54 i 57 ~ — analyze and dis- 
cuss the complications which arise when these forces are 
used to model structural relaxation in nanoscale struc- 
tures. It has been shown that such electrostatic forces are 
not conservative 5 -^— and may produce unstable or even 
run-away time evolution for the nuclei configuration . 48 ' 58 

The electrostatic forces can do work and will there- 
fore contribute to the total internal energy. However, the 
NEQ thermodynamic theory permits a quantification of 
the breakdown of the traditional HF theorem: 
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c/ LS 



Tr {(# 



H) 



Pls ■ 



(147) 
(148) 



The changes in the total internal energy content is not 
exclusively set by the work performed by electrostatic 
forces. 

Before proceeding, I stress that the lack of variational 
property for the total internal energy content is mirrored 
by a corresponding behavior for the Gibbs free energy 
state function. Taking the formal derivative, one obtains 



d 



, dY LS t 
Tr{(y LS 



Pls 



_d_f>-P(H-Yi,s) 



(149) 
(150) 



I note that for the widely used BTD method for 
electron-density calculations there does exist an argu- 
ments against exclusively implementing the electrostatic 
forces. Sec. VI. F shows that the BTD method for 
electron-density calculations includes a consistent de- 
scription of Gibbs free energy effects, namely by satis- 
fying the Friedcl sum rule in the BTD treatment of an 
effective SP scattering dynamics (even if this is not a 
complete treatment of MB interaction effects). This ob- 
servation suggests that it is natural to also include those 
Gibbs free energy effects when using the BTD result to 
implement forces on and relaxations of the classical nuclei 
under an assumption of adiabatic dynamics. 

For the equilibrium problem of surface-state mediated 
interactions between adatoms on Cu(lll), I note that 
the use of electrostatic forces (when evaluated in a GCE 
formal-scattering theory formulation 101 ) leads to a dis- 
agreement with experimental observations Using the 
GCE surface-state description, notation and concepts in- 
troduced in the subsection VIII. B, one readily obtains a 
formal GCE determination of electrostatic forces between 
the two adatoms 



-,es,GCE 
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dd 



dd\ , N 

e F aW 2 ,(e F ,<i;J F ).{152) 



Using formal Fourier-transformation theory, for example, 
as formulated by Lighthillj 1 ^ 5 - it is straightforward to ob- 
tain an analytical evaluation in the limit of asymptotic 
adatom separation, q F d > it. In particular, it follows 
that the asymptotical internal energy variation is entirely 
dominated byiSi 



AjV 2s (e,d>7r/g F ;c5 F ) 



-) K 

IT 



QFdy/y 

(ln(y)/7r + cot(<S F ) + i) 2 1 



[l + (ln(y)/7r + cot((5 F )) 
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(153) 

K=e/eF 



The leading-order GCE determination, (|151|) and 
(I153p . of electrostatic forces coincides exactly with the 
forces produced by the Gibbs free energy variation 
Y 2 s(d;5 F ) w e F AiV 2s (e F ; <5 F ). The asymptotic decay 1/d 
of either of the individual GCE interaction components 
differs from the 1/d 2 decay form derived in a previous 
perturbative canonical-ensemble calculation^ - and with 
experimental observations j 1 ^ The full GCE thermody- 
namics account of interactions and forces, next subsec- 
tion, is important even in equilibrium and even in the 
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limit of 1//3 — > 0, when using a (GCE) formal-scattering 
theory ; 67 ! 101 ' 102 ' 106 for example, in the analysis of surface- 
physics interaction problems.-^^ 



D. Grand-canonical ensemble thermodynamic 
forces 

I argue that for the standard open-boundary model- 
ing frameworks iof tunneling, it is more consistent and 
correct to describe (i.e., mimic from an electron QKA 
through use of a NEQ BOA) structure optimization for 
a working NEQ tunneling systems in terms of thermody- 
namic or effective forces 



LS,R; = 



d 



.0 L s(Mi;{Ri}), 



(154) 



here restated to emphasize the dependence on state vari- 
ables. 

In a thermodynamic description, the classical nuclei 
still respond to electrostatic forces but the adiabatic re- 
laxations (of very heavy classical nuclei) will also re- 
flect a tendency to seek (thermodynamically) favorable 
electron distributions. 67 An important part of these dis- 
tribution changes are specified by a correct evaluation 
of the Gibbs free energy effects (which exists even at 
(3^0) and which arc directly specified by the Friedel 
sum rule as formulated in the GCE thermodynamics de- 
scription. I stress that, (unlike in the equilibrium case,) 
the NEQ MB scattering causes actual changes in the 
electron distribution s 77 ' 82 ' 87 ' 126 ' 127 (thus significantly en- 
hancing the changes described by a direct variation in the 
SP potential) and that the extent and effects of such MB 
electron-distribution changes will also depend on the nu- 
clei positions*!! 3 - 

For the study of surface-state mediated indirect elec- 
tronic interactions between Cu(lll) adatoms , 101 ' 103 I 
note that the evaluation of GCE thermodynamics forces 
yields an agreement with experiments! - 3 - For the model 
description of the Cu(lll) surface-state interaction, and 
in the notation introduced in subsection VIII. B we have: 



jiGCE 



fi?s(d;fr) 



de'AN 2s (e';5 F ) 



(155) 
(156) 



-CF 



2sm(S F ) \ 2 sm{2g F d + 25 F 1 
7T / (q F d) 2 



The latter expression 10 - 1 applies in the limit of asymptot- 
ical adatom separations (q F d > tt) where it is possible 
to use the Lighthill analysis!^ for analytical evaluation 
of the asymptotic behavior. The nonpertubative result 
(|157l) agrees with the experimentally observed^ large-d 
interaction form. Also, the decay 1/d 2 of the oscillatory 
interaction is consistent with a perturbative canonical- 
ensemble evaluation!2°- I stress, however, that this is 



an interaction result which only emerges after correcting 
the internal-energy variation for the (much larger and 
much longer ranged) Gibbs free energy term (I144[) . In 
other words, it is the GCE thermodynamics forces (and 
not the GCE electrostatic forces) which in this model 
calculation 101 gives the best agreement with the experi- 
mentally observed adatom-interaction on Cu(lll)! - 3 - 

Finally, I observe that these thermodynamical forces 
are essentially already available at the level of the BTD 
method. In addition, there exist calculations which de- 
termines transport properties within ReNQSM also using 
more complete (but also more complex) representations 
of the MB interactions^! 



E. Adiabatic nature of thermodynamic forces 

The state function r2Ls(/ i i! {R-i}) is by construction a 
conservative potential for the forces (|154l) . For a given 
applied bias (set of external chemical potentials) we can 
also vary the positions to optimize the configuration with 
the lowest possible NEQ grand potential value. Use of 
the GHF theorem (|8ip yields the formal determination 
(fTT|) of the thermodynamic forces by the electron system 
on the nuclei. In contrast, there are no such arguments 
for the derivatives of the total-internal-energy or of the 
Gibbs- free-energy state functions, Eqs. ()148[) and (|150[) . 
when evaluated individually. 

In steady-state tunneling, the derivatives of the state 
functions are formally related by 
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S[p] 



Sp 



{n hS [p}-(U LS [p}-Y hs [p})}. (158) 



As discussed for the general time-dependent tunneling 
case in the previous section, the minimization of the ther- 
modynamic grand potential reflects a maximization of 
the entropy tempered by the correct evaluations of the 
NEQ value of the total-internal and Gibbs free energies. 

The fact that a GHF theorem (|5Tj) does holds for the 
thermodynamic forces (|154|) implies that they are, in fact, 
adiabatic. This follows because implementing the ther- 
modynamic forces (|154l) with a sufficiently slow (exter- 
nal) nuclei dynamics allows the GHF to automatically 
ensure that the system is always kept at the local state- 
function value Hls value. By the extremal nature of the 
thermodynamic grand potential it therefore also follows 
that the density matrix must then be the actual solution 
of the NEQ tunneling problem (in the configuration to 
which the system has evolved). 

The adiabatic relaxations by the thermodynamic 
forces will proceed and further maximize the entropy 
of the combined electron-nuclei systems until, eventu- 
ally, the system reaches a globally (or locally) optimized 
maximum-entropy morphology. At the same time, the 
foundation in the exact QKA ensures that the GCE ther- 
modynamics optimization is at always performed subject 
to an stringent adherence to conservation laws. This is 
exactly what we want from a nonempirical theory and 
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modeling of relaxations in steady-state tunneling sys- 
tems. 



generic tunneling systems will relax to a steady state 

] 3C ^ lx ^ or 25,28,28,67,74,76,77,93,95,96 



F. Conserving adiabatic forces and stability 

The demonstration of a conservative nature of adia- 
batic thermodynamic forces also has consequences for the 
consistency of making a steady-state transport assump- 
tion. The finding of a conservative nature of the adiabatic 
forces implies that one cannot (contrary to recent discus- 
sions) expect an abundance of possible current-induced 
instabilities, at least not without extending the study to 
also include an electron-vibrational coupling. 

Of course, there are systems that do not reach a steady 
stated 6 , and for which the thermodynamic state function 
and adiabatic forces cannot be defined. However, in spite 
of the ability of the nonconservative electrostatic forces 
to do work, I find that there is no intrinsic property of 
adiabatic NEQ forces (described in the formally exact 
electron QKA) which would necessarily drive a poten- 
tially steady-state system away from a stable morphol- 
ogy. Specifically, if we succeed in moving a nuclei in a 
tunneling system along a closed path (subject to the adi- 
abatic property which is suggested by the large difference 
in nuclei and electron masses) we can expect no accumu- 
lation or depletion to arise in any of the thermodynamic 
state functions. 

This robustness, the argument that a steady-state as- 
sumption is generally consistent, is a simple consequence 
of the conservative nature of thermodynamic forces; an 
argument is nevertheless included for completeness. The 
initial state (before the transformation around a loop) 
is uniquely characterized by the solution density matrix 
Pls and potential value Uls- Following any transforma- 
tion around a closed loop we obtain a new density matrix 
pf and a new thermodynamic potential value Cl(^ s . How- 
ever, with adiabatic transformations one can ensure that 
Flsr wm always characterize the potential changes and 
that 
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By the extremal property of the thermodynamic grand 
potential it therefore also follows that we have returned 
to the exact same density matrix, p' = /5ls- in other 
words, the conserving nature of the thermodynamic 
forces ensures that all state functions must attain the 
exact same value as the initial state functions, 
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The argument above applies for investigations of the 
electron-dynamics problem, outside Joule heating or 
other possible current-stimulation effects arising with 
an inclusion of an electron-vibrational coupling i 45 i 61 
I find that the conservative nature of steady-state 
thermodynamic forces provides a robustness which 
is consistent with the common expectation that 



G. Variational evaluations of adiabatic changes in 
thermodynamic state-function values 

I conclude the discussion of the here-defined conserving 
thermodynamic forces by completing a formal determina- 
tion of changes in the steady-state entropy, total-internal 
energy, and Gibbs free energy contents which results from 
an adiabatic relaxation. 

The key observation is that the operator analysis of 
Mermin, summarized in Appendix A, can be used to ob- 
tain the parametric derivatives of the unnormalized NEQ 
solution density matrix, exp[— j3(H — Ils)]- For example, 
use of (|A10|) and (|A13[) directly leads to an evaluation of 
the entropy content, 



as 
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-j3 2 ([(H-Y hs )-(U hs -Y hs )] 



d(H-Y LS ), 



(161) 

a form which is entirely specified by differences in internal 
energy and Gibbs free energy operators (and differences 
in their values). 

There is a corresponding formal simplification in eval- 
uation of the expectation value of the Gibbs weighting 
factor, 



d 



(Uts - Y LS ) 



5R"r Ls+ / Ls 



(162) 



Using a cyclic permutation of operators in the trace yields 



1 dS LS 
P 2 dR t 

and it follows that 
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It is also possible to complete a formal evaluation of the 
derivatives of total-internal energy and Gibbs free energy 
(also when evaluted in isolation) thanks to the nature^ 7 - of 
the Gibbs free energy or electron-redistribution operator. 
In steady state tunneling the Gibbs free energy must, of 
course, also be time independent. The time evolution is 
naturally formulated in terms of the Liouville operation 
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The Gibbs free energy operator commutes with the in- 
teracting steady-state Hamiltonian H = H(t = 0) in 
steady-state tunneling problems. As pointed out by 
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Hershfield, 67 the role of Yls is to redistribute the elec- 
trons among the MB eigenstates of H, Ref. [fj^- Using 
the commutator result, [H, Yls] = 0, together with Mer- 
min's observation, (|A13|) . leads to the formal evaluations 
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of the total internal energy and Gibbs free energy deriva- 
tives, respectively. 

Finally I note that the reformulation of the derivative 
of the expectation value of the Gibbs weighting factor 
(|164l) can also be used to obtain an alternative expression 
of the here-defined thermodynamic forces, 



-pGCE 
r LS,R; 



35- ~ ^ LS ^s + (/V2)(gf-(g - Yls) 
1 + /3(f/ LS - Y ls ) 



. (168) 

This alternative force expression suggests a relation be- 
tween the here-presented variational thermodynamic the- 
ory and the variance principle explored in Ref. l94l 



straightforward generalization of Mermin's previous re- 
sults for equilibrium thermodynamics^ The paper fur- 
thermore demonstrates that using the uniqueness of den- 
sity proof of the LSC-DFT— permits the formulation of 
the exact NEQ thermodynamic theory in terms of den- 
sity functionals. 

I derive exact thermodynamic forces, that is, effective 
forces which describe the behavior of the ensemble that 
characterizes the NEQ tunneling system. The thermo- 
dynamic forces are different from the electrostatic forces, 
satisfy a GHF theorem, and are thus explicitly conserva- 
tive. The paper demonstrates that these themodynamics 
forces are adiabatic. 

In conclusion, I propose that use of the here-derived 
conservative, adiabatic thermodynamic forces constitutes 
a more correct and natural implementation of the BOA 
for use with a NEQ electron QKA of tunneling. This 
because the GCE thermodynamics forces are consistent 
with the Friedel sum rule and because the adiabatic na- 
ture suggests that the thermodynamic forces will provide 
a better modeling of the actual relaxations of morphology 
which arise under operating condition. 
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IX. SUMMARY 

A complete NEQ theory of tunneling must be able to 
reliably describe the structure of a tunneling system un- 
der operating condition and, in particular, to determine 
transport-induced morphology changes. This has been 
a challenge because in actual NEQ QSM calculations, 
the focus is on the electron behavior. Previous investi- 
gations have used electrostatic forces in a description of 
such relaxations. For the standard (open-boundary and 
thus infinite) model framework of tunneling, I observe 
that the need to include effects of the Gibbs free energy 
variation and a a possibility for a complex entropic be- 
havior in fully interacting steady-state transport 67 makes 
it natural to move to a GCE thermodynamic description 
of current-induced, forces and relaxations. I find it im- 
portant to consider and evaluate an additional force con- 
tribution, one which arises from the sensitivity of the 
NEQ electron distribution to the scattering potential, 
i.e., the nuclei positions^ This force contribution de- 
scribe Gibbs free energy effects which are essential to 
satisfy the Friedel sum rule. 

To deliver a description of such effective forces, this pa- 
per recasts a rigorous GCE thermodynamic description 
which is based on formal (electron) QKAi£ and which 
is closely related to the exact steady-state ReNQSM. 67 
Specifically, this paper formulates a functional for the 
NEQ thermodynamic grand potential. In part, this re- 
sult is made possible by a definition of an operator for 
the NEQ Gibbs free energy under general NEQ tunneling 
conditions. The paper demonstrates a variational prop- 
erty of the NEQ thermodynamic grand potential as a 
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Appendix A: Mermin-type variational 
thermodynamics for nonequilibrium tunneling 

This section presents a summary of the proof of the 
variational nature of the NEQ thermodynamic grand po- 
tential and of the global-minimum property of the ex- 
act (solution) NEQ collision density matrix. The argu- 
ments are in many ways a straightforward generalization 
of Mermin's proofed for corresponding properties of the 
equilibrium thermodynamic theory. The generalizations 
applies for time-independent or time-dependent interact- 
ing tunneling systems under NEQ conditions defined as 
a collision problem, Sec. IV. The proof of the variational 
property is only summarized here with weight given to a 
detail that differs slightly from the equilibrium case. 

The first observation is that it is sufficient to proof the 
extremal property for any given time and hence I shall 
suppress all mentioning of subscripts 'col', 'sc', or 'LS' (as 
well as possible temporal argument) in the following. A 
more compressed formulation of the extremal properties 
of the NEQ thermodynamic grand potential, the NEQ 
variational property can thus be formulated 

n[p] > n[p \ = Q (Ai) 
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where (for a Hermitian Xq) 
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O = -^ln(Tr{e-^}), 
n[ P ] = Tr{/5[X +/?- 1 lnp]} 
and where trial NEQ density matrices 



P = 



Tr{e-/3*} 



(A2) 

(A3) 
(A4) 



(A5) 



are given by any Hermitian operator X. 

As in the equilibrium theory, it is advantageous to also 
define a set of intermediate density matrices 



X\ = Ao + AA, 



(A6) 
(A7) 



where A = — (1//3) ln(/5) — Xq. This family of possible 
NEQ density matrices (|A6I) form a line between po and p 
where each member formally represents (in the notation 
introduced in Sec. VI) a solution of the thermodynamic 
grand potential 



Tr{e-^} 
-iln(Tr{e-^})^n[p A ]. 



(A8) 



I note that, unlike in the equilibrium case, a sepa- 
rate argument (for example, V-representability) would 
be needed to conclude that any given member, p\ (of the 
line of density-matrices) corresponds to an actual den- 
sity matrix solution for a physical NEQ tunneling sys- 
tem. However, such a representability assumption is not 
of any consequence for the proof and discussion of the ex- 
tremal nature of the NEQ thermodynamic theory. The 
essential observation is that all of the intermediate NEQ 
density matrices (|A6j) are still specified in the form given 
in Eq. (|A3|) . 

The central step in the proof and argument is a demon- 
stration of the variational property, 
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that is, a compact formulation of the generalized 
Hellmann-Feynman theorem (|81l) . Here the key obser- 
vation is the general operator identity^ 
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which, of course, applies irrespectively of the nature of 
the statistical ensemble, for NEQ and equilibirum condi- 
tions alike. Introducing also (0)\ = Tr{p\0}, one can 
thus formulate the parametric density matrix changes, 
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Mermin also makes the observation that a cyclic per- 
mutation of operators permits a significant simplification 
in the evaluation of expectation values, for example^ 



<A A (/3')>x = <A)x. 



(A13) 



Inserting this result in (|A12[) completes an explicit 
demonstration of the generalized Hellmann-Feynman re- 
sult (|A9]). 

Finally, from (IA10[) also follows a proof for the ex- 
tremal nature of the thermodynamical grand potential. 
It proceeds by the reformulation of the potential deriva- 
tive 



(A14) 



= A / d/3'|(A A (/3/2)-(A) A )| 2 . (A15) 
Jo 

As in equilibrium thermodynamics, this results is suffi- 
cient to ensure that the difference fl[p] — 0[po] is positive 
and cannot vanish except when p = po- 

In summary, the operator analysis of Mermin 16 is suffi- 
cient to ensure (not only for the equilibrium case but also 
for the present NEQ thermodynamics) that the thermo- 
dynamic grand potential is both variational and extremal 
at the exact solution /5q. 



Appendix B: Linking LSC-DFT and the quantum 
kinetics accounts of tunneling 

The choice of the continuum Caroli partition scheme 
is designed to allow a combination of the LSC-DFT 28 
with the present general NEQ thermodynamic account. 
To accomplish that, it is necessary to also ensure that 
the choice of details in the adiabatic turn (for the elec- 
tron QKA) is made such that it permits the LSC-DFT 
proof of uniqueness of density. The underlying issue is 
purely technical, but for completeness, I include here a 
mathematical analysis. 

While Ref. [28| never explicitly stated this, the for- 
mal LSC-DFT proof of uniqueness of density does, in 
principle, require that the adiabatic-turn-on factor a r) (f) 
vanishes not only at t — > — oo, but also for all times 
t < to where to — > — oo. In a strict mathematical sense, 
this is necessary because the LSC-DFT proof implicitly 
assumes one can work with a Taylor expansion in the 
temporal evolution. With the simple adiabatic turn on 
[an(t) = exp(?7i)] that we normally employ in the QKA, it 
follows that all Taylor-expansion coefficients of the tem- 
poral evolution vanish in the limit t — ¥ —oo. It is clear, 
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that alone an extremely soft modification of the adiabatic 
turn on is required to ensure a rigorous mathematical 
description in the LSC-DFT analysis. It is also clear 
(because of the entropy flow) that such an extremely 
soft modification of the adiabatic turn on can have no 
consences for the tunneling dynamics at relevant time, 
for example, as expressed in the steady-state LS solution 
to the MB collision problem. The argument can be made 
formal. 

To be specific I investigate adiabatic turn off which [in- 
stead of the standard LS or QKA form, a v (t) = exp(r/t)] 
is identically zero until some far away time to —> — oo. I 
assume that a v (t > to) is real, monotonously increasing, 
and satisfies < a v (t) < 1, with a bounded derivative 
(darj/dt) < r]ari(t). An explicit choice of slightly modi- 
fied adiabatic turn on factors can be expressed 



a n {t) = exp(rjt)b tl (t;t2,t^) 



(Bl) 



where I have introduced a monotonously increasing pos- 
itive function *q, t^) with b v (t < t^) = 0, b„(t > 
tl) = 1 and db v {t)/dt < 77. I choose t\ = 
(with a — > 00) and t r ' = t^/r/ 2 . I note that this set of 
choices removes any possible t — —00 singularities in the 
LCS-DFT proof of uniqueness of density and I proceed 
to demonstrate that this modified adiabatic turn on still 
leaves all details of the formal LS solutions in the QKA 
unchanged (in the relevant limit 77 — > + ). 

The important fact is that da^fy/dt remains essen- 
tially specified by ricxp(r]t) since this factor, in turn, de- 
termines the denominator in the formal MB LS solution 
(jltl)) . I introduce 



m = 



-i(E i -H+iij)t 



dt'Ut') 



i(H-E ( +ir])t 



E e -H 



(B2) 



(B3) 



and note that the boundary conditions, i.e., the original 
choice of adiabatic turn on exp(?yi), enters the MB LS 
solution through the formal evaluatio n 109 ' 147 



I* 



J +) )=7^(0)|$ C 



(B4) 



Using instead the present collision picture with slightly 
modified adiabatic turn on (|B1() , the boundary conditions 
instead enters the formal MB solution through 



f dt> 7e(OMO = hi* = 0) 

J —00 



Bdv) (B5) 



t? 



B.UD = I ,,,'^(1/)^), (B6) 



as obtained by partial integration. The resulting LS so- 
lution therefore has the form 



|^(+)) = v At(ri)\$t). _ (B8) 

However, the relative normalization of the additional 
term in the modified MB LS solution (IB8|) is strongly 
bonded, 



< 77 / 1 dt' 
ItV, 



< Bin) |*[ +) ) 



^'"'^(0)1^) , (B9) 
(BIO) 



by the (state-independent) evaluation 



B(rj) = exp«) = exp[-(l/r y ) a ] -► 



(Bll) 



The choice t'{ = ?^( Q + 1 ) with a — > 00 certainly ensures 
that B{rj) and hence the modification (|B8I) of the LS so- 
lution become irrelevant in the limit 77 — > + . 

In summary, I find that the LSC-DFT analysis 2 - can be 
adapted to provide an independent proof for uniqueness 
of density while leaving the formal LS solutions of the 
MB collision problems unchanged. This opens for a den- 
sity functional formulation of the NEQ thermodynamic 
quantities in a description which resembles the equilib- 
rium case presented in Ref. [l6l 
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